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Summary  of  the  Effort 

The  research  effort  emphasized  system  identification,  decentralized  control  design, 
and  their  applications  to  smart  structures  and  flexible  pointing  systems.  The  research 
effort  had  both  theoretical  and  experimental  ingredients.  On  the  theoretical  side,  new 
system  identification  techniques  for  nonrational  systems,  which  may  be  either  min¬ 
imum  phase  or  nonminimum  phase,  were  developed.  The  identification  techniques 
complemented  our  efforts  on  theoretical  modeling  of  systems  with  embedded  and/or 
surface  mounted  smart  materials.  Having  attained  a  theoretical  model  for  systems 
at  hand,  robust  linear  and  nonlinear  control  techniques  are  utilized  for  performance 
enhancement.  The  control  strategies  are  robust  to  structured  and  unstructured  per¬ 
turbations  of  the  system  dynamics.  The  control  design  methodology  are  based  on 
a  decentralized  framework.  This  provides  local  decision  making  through  the  smart 
structures  when  multiple  sensors  and  actuators  are  present.  On  the  experimental  side, 
we  fabricated  several  smart  composite  tubes  with  piezoelectric  sensors  and  actuators. 
One  such  tube  (with  surface  mounted  piezoceramic  sensors  and  actuators)  was  used 
to  replace  part  of  the  existing  weapon  pointing  setup,  namely  ATB-1000,  at  the  Army 
Research  Laboratory  at  Picatinny  Arsenal.  Furthermore,  we  have  developed  similar 
experimental  setups  including  a  two-link  flexible  arm  with  piezoelectric  sensors  and 
actuators  at  the  Control/Robotics  Research  Laboratory  (CRRL)  at  Polytechnic.  The 
organization  of  the  final  report  is  as  follows:  a)  the  decentralized  control  design  frame¬ 
work  developed  is  reported  first,  b)  the  system  identification  for  non-minimum  phase 
systems  is  presented  next,  c)  lastly  a  brief  outline  of  our  experimental  effort  is  given. 

Keywords:  Nonlinear  and  Decentralized  Control,  Adaptive  Control,  System  Identi¬ 
fication,  Non-minimum  Phase  Systems,  Smart  Structures. 
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1.  Decentralized  Adaptive  Output  Feedback  Design  for  Large-Scale 

Nonlinear  Systems 

1.1  Introduction 

Recently  adaptive  control  techniques  have  been  utilized  for  decentralized  control 
of  systems  with  uncertain  interconnections  bounded  by  first  order  polynomials [1-7]. 
The  first  available  result  on  decentralized  adaptive  state  feedback  control  for  higher 
order  interconnections  was  reported  in  [8].  A  major  structural  restriction  in  [7,8]  is 
that  to  ensure  global  stability,  the  interconnections  must  lie  in  the  span  of  the  input. 
Results  on  decentralized  adaptive  control  for  linear  systems  with  subsystems  having 
arbitrary  relative  degrees  are  presented  in  [1-4].  For  higher  order  interconnections, 
the  strict  matching  restriction  was  relaxed  in  [9],  where  decentralized  adaptive  state 
feedback  controllers  were  designed  using  tools  from  centralized  nonlinear  control  [10]. 

Relatively  few  results  on  decentralized  output  feedback  are  available  in  the  liter¬ 
ature,  and  still  fewer  for  the  adaptive  case  [5-7].  Recently,  based  on  adaptive  output 
feedback  design  in  the  centralized  case  [11],  a  decentralized  adaptive  regulator  was 
designed  for  in  [3]  for  systems  of  the  observer  canonical  form.  In  all  these  output 
feedback  schemes,  the  interconnections  are  assumed  to  be  bounded  by  Lipschitz  type 
bounds. 

In  this  report,  global  decentralized  adaptive  output  feedback  controllers  are  de¬ 
signed  for  large  scale  nonlinear  systems  of  the  output  feedback  canonical  form.  We 
consider  both  parametric  and  dynamic  uncertainties  in  the  interconnections  under 
the  assumption  that  the  interconnections  are  bounded  by  unknown  pth  order  poly¬ 
nomials.  In  addition,  we  also  consider  bounded,  unmeasurable  disturbances  entering 
the  system.  The  design  starts  with  the  construction  of  decentralized  filters  for  state 
estimation  for  each  subsystem.  The  estimation  error,  however,  is  not  exponentially 
convergent  due  to  interconnections  and  disturbance  terms.  Utilizing  techniques  from 
[10,12],  a  systematic  procedure  is  developed  for  design  of  an  adaptive,  decentralized 
controller  which,  through  injection  of  nonlinear  damping^,  guarantees  global  uniform 
boundedness  of  the  tracking  error  and  all  the  states  of  the  closed-loop  system  to  a 
compact  set,  and  also  provides  good  disturbance  rejection  properties.  The  size  of  the 
compact  set  can  be  made  arbitrarily  small  through  proper  choice  of  the  control  gains. 
A  higher  order  composite  Lyapunov  function  is  employed  to  prove  global  decentral¬ 
ized  stabilization.  Furthermore,  while  the  majority  of  output  feedback  schemes  in  the 
centralized  case  can  handle  only  linear  parametric  uncertainties,  the  proposed  scheme 
incorporates  both  parametric  and  dynamic  uncertainties  in  the  interconnections,  and 
the  uncertain  parameters  may  appear  nonlinearly  in  the  vector  fields.  In  addition,  the 

^The  nonlinear  damping  terms  designed  in  this  report  are  different  from  those  obtained  by  Kanel- 
lakopoulos  et.  al  [10]. 
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proposed  scheme  requires  adaptation  of  oniy  one  scalar  parameter  for  each  channel 
in  contrast  to  [11],  where  the  number  of  adaptations  is  the  same  as  the  number  of 
uncertain  parameters.  Furthermore,  no  controller  redesign  is  required  if  additional 
subsystems  are  brought  on-line  or  taken  off-line.  Global  asymptotic  properties  are 
attained  for  the  case  where  the  objective  is  regulation. 

1.2  The  Class  of  Large-Scale  Nonlinear  Systems 

We  consider  a  large-scale  continuous  time  affine  nonlinear  system  comprised  of  N 
channels  (agents)  given  by 

N 

X  =  f{x) +  Yj9i{^)'^i  + 

i  =  l 

Vi  =  ^i(^)  l<i<^  (1-1) 

where  x  G  7^^  is  the  state  vector  for  the  overall  system,  Ui  G  and  yi  G  'R  are 
the  input  and  output  for  the  ith  subsystem,  and  uj{t)  is  a  bounded  unmeasurable 
disturbance  signal.  The  vector  fields  /,  p,  gi^  and  the  function  hj  are  assumed  to  be 
smooth  on  with  /(O)  =  0,  p(0)  =  0,  hj{0)  =  0,  and  pi(0)  ^  0.  Only  the  output 
yi  is  available  for  feedback  to  the  zth  control  channel. 

The  decentralized  scheme  proposed  in  this  report  is  applicable  to  a  class  of  large- 
scale  nonlinear  systems  which  are  transformable  via  a  global  diffeomorphism 

^  =  [211  .  .  .  2i«i2i2  ...Zm...  $(0)  =  0 

into  the  output-feedback  canonical  form: 

Zii  =  Zi2  +  (fiiivi  ...yN)  +  CniVi  ■  ■  ■  yN)^{t) 

^i,pi  —  ^i,p,+i  d"  'r’itPiiyi  •  •  -  yN^  d"  G',Pi(yi  •  •  •yN)^{t)  bi,Ki-pi^i{yi)ui 

=  y^i,Ki{y\---yN)  +  ii,f^Xy\---yN)‘^{t)  +  bio8i{yi)ui 

yi  =  Zii,  1  <  i  <  A''  (1.2) 

where  1  <  z  <  are  the  observability  indices  with  ^nd 

Ylf=\  Ki  =  Af,  and  . . .  b^o  are  unknown  constant  parameters.  The  coordinate- 

free  necessary  and  conditions  for  transformation  of  (1.1)  to  (1.2)  are  given  in  [14]. 

We  consider  both  parametric  and  dynamic  uncertainties  in  the  terms  (pij{yi  . . .  j/iv) 
and  ^ij{yi  . .  .J/at),  Specifically,  we  assume  that  are 

unknown  except  that  they  are  bounded  by  an  unknown  Pijth  order  polynomial  in 
outputs,  i.e., 
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Vi  3  N  Pij  N 

<  EEfSibfi  <  EE4illrf  (i-s) 

k-O  l=l  k=0  1=1 

with  (!iji  and  i^^ji  unknown. 

Remark  2.1:  Note  that  (fij  and  may  be  nonlinearly  parameterized  in  the  un¬ 
known  parameters  provided  condition  (1.3)  is  satisfied.  Also,  only  the  quantity 
p  =  maxi<KjV;  i<j<«i  Pij  needs  to  be  known. 

To  avoid  losing  insight  of  the  underlying  idea  of  this  report,  we  assume  that  the 
parameters  bij,  0  <  j  <«,  —  /),,  1  <  i  are  known.  The  case  of  unknown  bij  can 
be  treated  identically  as  in  [11]  where  the  following  assumptions  are  made: 
Assumption  (2.1):  The  sign  of  bi^^i-pi  is  known. 

Assumption  (2.2):  The  polynomial  Bi{s)  =  +  •  •  •  +  bus  +  6,o  is  strict 

Hurwitz  (exponential  stability  of  the  zero  dynamics). 

Assumption  (2.3):  <5j(t/,)  ^  0, 

Assumption  (2.4):  The  reference  signal  yi^rej  and  its  first  pi  derivatives  are  known 
and  bounded  with  piecewise  continuous. 

For  known  bij,  we  assume  assumptions  (2.2),  (2.3),  and  (2.4)  to  be  true. 

1.3  Decentralized  Adaptive  Design 

The  decentralized  adaptive  control  scheme  is  designed  for  the  transformed  system 
(1.2).  Since  only  the  ith  output  (yi)  can  be  utilized  for  feedback  to  ith  control 
channel  (u,),  we  need  to  design  an  observer  to  estimate  the  remaining  states  of  the 
ith  subsystem.  However,  an  observer  with  linearizable  error  dynamics  cannot  be 
designed  in  this  case  due  to  parametric  and  dynamic  uncertainties  in  the  (fij  terms 
in  (1.2).  Even  if  the  exact  form  of  ifij  is  known,  in  the  decentralized  context,  the 
presence  of  unmeasurable  outputs  from  other  subsystems  yk,  I  k  <  N;k  ^  i  makes 
it  impossible  to  obtain  a  linearized  error  dynamics.  However,  as  will  be  shown  in 
the  ensuing  design  procedure,  the  observer  error  can  still  be  made  to  converge  to 
zero  asymptotically  using  nonlinear  damping  terms  in  the  control.  The  systematic 
design  procedure  presented  here  employs  tools  from  the  nonlinear  design  toolkit  by 
Kanellakopoulos  et  al.[10]. 

The  ith  subsystem  (1.2)  can  be  written  as; 

=  AiZi  -h  kiyi  +  (pi{yi  ...yi^)  +  (i{y-i  . . .  ?/v)<^(t)  +  bi8i{yi)ui  (1.4) 

where 

—  :  I  ^  ki  —  ;  •)  —  0  ‘  *  0  bi^tzi—pi  *  *  *  biQ  , 

0  *  ■  *  0 

X*  "  lx' 

J  ,  —  ^il  •  •  • 


The  gains  ki  are  chosen  such  that  Ai  is  a  strict  Hurwitz  matrix. 

We  define  the  estimate  of  the  states  Zi  as 

Zi  =  AiZi  +  kiVi  +  bi6i{yi)ui.  (1.5) 

The  observer  error  Ci  —  Zi  —  Zi  is  given  by 

ki  =  Aie  +  . . .  ?/;v)  +  ^i{yi  ■  ■  ■  (1-6) 


Since  Ai  is  a  strict  Hurwitz  matrix,  given  a  positive  definite  matrix  Qi,  there  exists  a 
positive  definite  solution  Pi  to  the  Lyapunov  equation 

AjPi  +  PAi  =  -Qi.  (1.7) 


In  the  following,  we  assume  uniform  relative  degree  for  each  subsystem  {pi  =  p, 
1  <  I  <  iV)  in  order  to  reduce  notational  complexity.  The  case  where  the  relative 
degrees  are  different  for  each  subsystem  can  be  treated  similarly. 

The  ensuing  systematic  design  procedure  involves  at  each  step,  augmenting  an  inte¬ 
grator  and  the  design  of  a  virtual  control  law  to  stabilize  the  augmented  system,  until 
the  actual  control  appears  at  the  pth  step. 

Step  0:  Define  the  tracking  error  for  the  fth  subsystem  as  Xi\  —  “  Vi, ref  =  Vi  ~ 

Tio{yi,ref)-  Expressing  z,2  in  terms  of  its  estimate  as  Zi2  =  2,2  +  ^i2i  we  obtain 
Xil  =  ^i2  +  ‘fil{yi  ■■■yN)  +Cil{yi  ■  •  •yw)w(t)  -  in, ref  +  e,'2,  1  <  f  <  iV.  (1.8) 

The  following  virtual  decentralized  control  law  2, 2  for  the  fth  subsystem  (1.8)  is  pro¬ 
posed 

Zi2  =  -TTiiXil  -  0iXil  -  neiixil  +  xlr^)  +  in, ref  =  Tii(Xil,  A)  +  in, ref  (1-9) 


where  Qi  >  0,  (3i  is  a  time  varying  adaptation  gain  to  counter  the  effects  of  the 
interconnections,  tvh  is  a  constant  gain  whose  magnitude  will  determine  the  ultimate 
bound  on  the  tracking  error,  and  7^  is  a  constant  gain  given  by 


li  =  (\  +  Np„)t,e  =  +  2Np.Mp  +  \)(2p  +  1) 


(1.10) 


Since  2,2  is  not  the  actual  control,  define: 

Xi2  “  ^i2  ’rji(Xil)^i)  yi,ref‘  (EH) 


Therefore,  Xii  =  Xi-2  +  r,i(xii,  A)  +  <Pii{yi  . . .  2/w)  +  6i(?/i  •  •  •  yN)^{t)  +  ^i2- 

Let  A*  he  the  desired  value  of  the  control  gain  A  to  counteract  the  effects  of  the 

interconnections.  The  following  adaptation  for  A  is  used 

A'  =  r,^if:Axn -^lUA  (1.12) 

k=l 
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where  ai  >  0  incorporates  the  standard  “cr  modification”  used  in  robust  adaptive 
control  to  avoid  parameter  drift.  The  true  value  of  /?*  required  to  counter  the  effects 
of  the  interconnections  is  obtained  in  the  final  step. 

Consider  the  following  composite  Lyapunov  function  for  the  sub-subsystems  (1.8) 


k=l 


(1.13) 


where  Fj-  is  a  positive  constant.  Differentiating  Vo  along  the  trajectories  of  Xii,  we 
obtain 


^0  =  -ejQiti  +  2eJPi  {ipi{yi . . .  y;v)  +  (i{yi  ■  •  •  J/)v)w(0)  +  ^[Xi2 

i=l  I  k=l 

+'Til{Xily  Pi)  +  ^i2  +  ■  •  •  Vn)  +  CiliVl  ■  ■  ■  yN)<^{t)]  +  ‘^^7^  {Pi  ~  Pi)Pi^ 

<  Y1  {-^mi«((5i)||etj|^  +  2||eiH||PiH  . . .  ?/n)||  +  UiiVi  ■  •  •  2/w)c^(0ll) 

2  =  1 

—  XI  xll  +  X  ‘^^SiXil  ^[xi2  —  7«Tt(Xjl  +  Xii  ~  PiXil  +  ^12] 

k—l  k=:l 

+  ■■■  !/«)ll  +  II6'i(!/i  ■  ■  ■»)"(<)ll) 


+2(ft-A-){rr'/3,-t.,y;tx?.‘}  . 

k=i  ) 

Utilizing  the  bounds  (1.3)  and  inequalities 
2ab  <  +  6^; 


Ykz=\ 

^  2 


\k=l 


\k=l 


Xkl  <2'=X 


\k 

Cli\  , 


(1.14) 

(1.15) 

(1.16) 

(1.17) 


.1  =  1 


i  =  l 


the  terms  in  (1.14)  can  be  bounded  as 


N  V 


N  P 


IIw.(!/i...w)II<EE24,IIx.iII“-' 

2=1  k=l 


Pti  N 

E  Ecfellwf- 


A:i=l  j=l 


(Using  (1.3)) 


N  Pil  N  (  p 

XXX2  X^^*llxair'=-^ 

i=lfcl=lj  =  l  IA;=1 
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N  pa  N  [/  P  V  ,  .2  1 

^  E  £  E +(cy  iiyiir''^  [  (Using (1.15)) 

i=l  ki=l  j=l  [  \k=l  /  ) 

N  Pi\  N  p  p  N  Pii 

<  E  E  E{E*:"??}{E +  EE'in*(ll».«/ll  +  ll»illf‘ 

2=:lA:i=li=l  k=l  k=l  i=l  k=l 

(Using  (1.16));  <ii,i  =  E  fe)' 


i=i 


N 


P  N  Pii  N  p 

2=1  k=l  2=1  k=l  2=1  k=l 


(Using  (1.17));  ^  N pn  Y,  = 

k~l 


Npi\p{p  +  l)(2p  +  1) 


Similarly, 


(1.18) 


EE2tPillx..ll“-'Ki(!/i..-»Mf)li 

2  =  1  ^=l 

N  p  N  pii  N  p 

i=l  A:=l  i=l  A:=l  2=1  k=l 


\2k 


N  2 

with  7ii  as  in  (1.18);  vuk  =  X] 

J-1 


(1.19) 


E E <  E2 (f: i?.||Xiill“-'|  l|£,ll 

2=1  k=l  2=1  l/:=l  ) 

^  +  l)(2p  +  1)  IkilP,  (1.20) 

2=1  A;=l  2=1 


2  =  1 
|2 


E2||Eill  (l|P.llll9.(!/.  . . .  »)ll)  <  E  l|£,f  +  E  l|Pini?>,(!/i . . .  w)IP 

2=1  2=1 

<  Eii£iir+EEiiPiriiv.s(!/i--.x»)iP 

7=1  2  =  1  j  =  l 

N  N  Ki  Pij  N  ^ 

<  Eli'.ir  +  EEEEll^’.ffe)  fell""  (Using(i.i7)) 

2  =  1  j  =  l  k-1  1-1 


2  =  1 
N 


< 


2  =  1 


l\  N  p  /V 

X  IPdP  +  X  X  d^k\\y^f’"]  where  d,k  =  X  ll^'lP  X  (&) 

i  =  l  k  =  l 


N 

C 

1=1 


J=l 
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(1.21) 


<  E I|C'II“ + E E +  EE 

j=l  i=l  k=l  i=l  k=l 


E21l«i|l(ll«IIK(!/.---w)"(OII)<EI|£.lP  +  EE2“«il!/i,wll“ 

izzl  i=l  t=l  k=^l 

+  X]  X]  fJik  =  ll-P/ir  'yi  (^maxl^Ui)  ■ 

i=l  k=l  1=1  j=l 


(1.22) 


With  adaptation  (1.12)  and  utilizing  the  bounds  (1.18)-(1.22),  and  noting  that  7^  = 
in  +  I  Efc=i  Vo  may  be  written  as 


l>o  <  E  I  -{A™i„(a)  -  3}l|e,||^  -  27r,i^,  E  xli  +  E 

i=l 


+‘^iiQi  E 


\\4k-2 


k=l  k—1 

P 


"^HqKxh  +  Xa  xli  ^  “  2cri(/5^  -  A*) A 


/c=l 


+  ^  +  ^tAr  +  diik  +  ^lz/c)(||xilir^  +  ||2/t,re/|l^^)  Qi  |  (1*23) 

k=l  k=l  J 


N 


<  E  -{A«.n(«.)  -  3)|leif  -  2^„e.  E  x“  +  E  2fci>.Xa‘'x.2  -  t'Kft  -  ft*)" 

i=l  I  /:=!  k=l 

+  E  2^^(difc  +  Vik  +  d-iik  +  niiA:)(||Xil|P*^  +  ||yi,re/||^^)  ~  Qi  E  ^Xf^'^) 


Inequality  (1.24)  is  obtained  from  (1.23)  by  completing  squares  of  the  {i^i  —  13*)  terms 
and  using  the  fact  that  the  resultant  of  the  first  two  terms  in  the  second  line  of  (1.23) 
is  negative. 

The  state  equations  for  xn-i  Xi2  given  by 


Xii 

Xi2 


=  Xi2  +  Tii(xii,  A)  +  v’iiiyi  •  •  •  yw)  +  (ii{yi  ■  ■  ■  yN)<^{t)  +  £*-2 


+  ki2{yi  —  Zil) 

5t,1 


dTi 


!l 


Qt 

{Xi2  +  ^^(Xtl,  A)}  -  ~ 


dx 


il 


dxil  . . . .  dl3, 

Wiiiyi  •  •  •  yyv)  +  61  (yi  •  ■■yN)oj{t)  +  £,2} 


—  2)3  T  ‘512(2/1’)  ^ili  ^i2i  (3ii  yi,refi  yi,Tef  i  yi,Tef) 

+^i2{xii,  ^i){(-i2  +  <fii{yi  •  •  •  yw)  +  ({yi  •  •  •  yN)^^{t)}- 


(1.25) 


Remark  3.1:  The  terms  {duk  +  n,u)|lxii|P  in  Vq  arise  due  to  the  interconnection 
and  disturbance  terms  appearing  in  the  state  equations  of  the  plant  (see  e.g.,  (1.8), 
(1.25)).  However,  as  can  be  seen  from  (1.24),  these  terms  are  in  the  span  of  the  control 
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and  can  be  countered  using  sufficiently  high  control  gain  j^*.  In  the  subsequent 
steps,  additional  terms  of  the  type  (OlIXulP^  will  appear  due  to  disturbance  and 
interconnections.  The  desired  value  of  (3*  to  counter  these  effects  is  assigned  in  the 
final  step.  Similarly,  possibly  destabilizing  disturbance  and  interconnection  terms 
arise  in  the  observer  error  dynamics  (1.6).  These  appear  as  [dik  +  Utfc)||Xji|P* 
the  Lyapunov  function  and  can  again  be  countered  by  injecting  nonlinear  damping 
through  /3*. 

The  design  procedure  can  now  be  made  recursive  as  shown  in  the  following  step: 
Step  m  (  I  <  m  <  p  —  2) :  \n  steps  0  to  m  —  1,  the  virtual  controls  th  . . .  Ti^rn  were 
designed.  In  the  m  —  1th  step,  Xi,m+i  is  defined  as  the  difference  of  5i,rn+i  and  Ti^m, 
i.e., 

X*,m+1  “  l3i^yi,rej  •  •  ■  Vi^ref^i  (1.26) 

The  dynamics  for  Xt,m+i  are  given  by 


Xi,m+1  —  ^i,m+2  T  ?t,m+l(yt?  ^il  ■  ■  •  Vi, ref  ■  ■  ■  Vi^ref  ) 

T^i,m+1  (j/i)  ■^il  •  •  •  ^i,m+l  •>  ^ii  Vi^ref  •  •  •  Vi^ref  )  {^i2  T  ■  •  •  J/iv) 

+^ii(yi...yN)a^(0} 


(1.27) 


where  the  notation  in  (1.27)  is  similar  to  that  used  in  (1.25).  We  consider  Zi^m+2 
as  the  virtual  control  for  the  (xn  ■  ■  ■  Xi,m+i)  sub-subsystems  and  use  the  following 
composite  Lyapunov  function  for  this  subsystem: 

K.  =  K,-. (1-28) 

l=\ 

Differentiating  along  the  trajectories  of  (xa  •  •  •  we  obtain 


N 

K.  <  E 

Z  =  1 


-{KiniQi)  -{m  +  2)}||e;f  -  ai0i  -  -  2TViiQi  Y,  xli  - 

k=l  j=2 

,m+l  •)  Vi, ref  •  •  •  yi,ref  )} 

“i"2Xt,m-t-l  (2/t  7  A7  Vi, ref  •  •  •  Vi^reJ  ) 

{^2  +  . . .  t/iv)  +  6i(j/i  •  •  •  2/n)^(0} 

+  X/  +  Vik  +  'fT^idlik  +  '^lik)}  {\\Xil\\‘^^  +  II  ||^^} 


k=l 


+cTi/3f  -213-g^^kx 


k=i 


(1.29) 


Note  that  for  m  =  1,  instead  of  the  term  2xi2Xti  ia  (1-29),  we  have  2x^2  Z)fc=i  kpiX^Y^. 
Utilizing  bounds  (1.3)  and  inequalities  (1.15)-(1.17),  the  terms  in  (1.29)  may  be 
bounded  as 
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N  N 

Y^2Xi,m+l^i,m+l^i2  <  +  ll^!'ll^) 

i=l  i=l 


(1.30) 


N 


N 


N  N  Pii 

'^Xi,rn+l'^i,7n’^l^il{yi  •  •  •  Vn)  ^  Xz,m+l^i,m+l  ”l~  X^  X^  (*^Ul)  WV^W 

2=1  j=l  k=l 

n  p 


2k 


2  =  1  2  =  1 

N  N  Pii 

~  XJ  X2,m+1^2,m+l 

2=1  2=1  k=l  2=1  /i;=l 


(1.31) 


N 


^  y  2Xt,m+l^i.m+l‘atl(i/l  •  •  •?/iv)^(0 


i=l 


n  p 


N  N  Pii 

^  2^^^^1i*i|bi,r-e/ir^  +  XlIZ  2^''t>iifc||x,l  | 

t=l  k=l  i=l  k=l 


2k 


(1.32) 


2=1 


where  duk^  vuk  are  as  in  (1.18)  and  (1.19).  Using  (1.30)-(1.32),  Vm  is  written  as 

p 


N 

Kn  <  E 

2  =  1 


-{XrmniQi)  -  (m  +  3)}|le,||2  -  <70,  ”  A*)"  "  27ra£>,  E  Xa  -  2 

k=\  3=2 

+2x  2, m+1  {-^2,771+2  "f"  Xim  i  •>  ^il  •  •  •  Vi, ref  •  •  •  Vi^ref  ) 

^ Xi^m+l'^i^rn+l  iVi^  ^il  ’  •  •  ^2,m+l7  A'?  Vi, ref  •  •  •  Vi^ref  )} 

+  ^2  +  Vik  +  (m  +  l){diik  +  ^l2^0}{|lx^l|^^  +  Il2/i,7’e/ii^^} 

k=l 

+cr,l3f  -  2/?-Si  E  txn  .  (1-33) 


k=l 


With  Zi^rn+2  as  the  virtual  control  for  the  (xn  •  •  •  Xi,m+i)  subsystem,  we  choose 
-2^2,777+2  ~  [^2,m+l X«,m+1  H”  X^m  ~t~  Q,m+l(j/2  7  -^21  •  •  •  A  7  yi,ref  •  •  •  yl^ref 


(m+1) 


“I" X2,m+l^i^7n-|-l  {yi^  >  A?  yi,ref  •  •  •  yi^ref 


A 


(1.34) 


”^t,m+l(yt5  -^il  •  •  •  -St.m+li  ^ii  Vi, ref  ■  ■  ■  Vi, ref  ) 

with  TTi^rn+l  >  0.  Define:  Xi,m+2  =  2t,m+2  -  Ti^m+liVi,  A,  Pi, ref  ■  ■  ■  2/1, Te/^*)- 

Then,  Vm  is  given  by 


N 


14.  <  E 


2  =  1 


m+1 


-{0in(Qi)  -  (m  +  3)}||ei|p  -  (7i(^i  -  -  2'KixQi  E  xfi  -  2  E 


A:=l 


i=2 


+  2X2,m  +  lX7,m+2  +  X^  2^^{dik  +  '^27:  +  (^  +  l)(^^l2A:  +  ^liA:)}  {  ||Xtl  11^^  +  ||2/2',re/ ||^^} 


k=l 
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(1.35) 


+cr,fif  -2l}’eiJ2l^Xa  ■ 

k=l 

Step  /9  —  1;  In  step  p  —  2,  Zi^p  was  derived  as  the  virtual  control  input 

Zii . . .  y^^ref  . . .  2/lrI/^))  foi'  the  {xn  . .  •  subsystem.  Define 

Xi,p  ~  ^i,p  ~  '^i,p-l{yi')  ^il  *  •  •  Viyref  •  *  •  Vi^ref  )* 

The  actual  control  input  Ui  now  appears  in  the  state  equation  for  Xi,p^ 

Xi,p  bi,K^~p^i{yi)'^i  T  ^i,p-\-l  T  ^il  •  •  •  ^i.p-)  A?  Vi, ref  •  •  •  Vi^ref) 

+w^^p{yi,  z,i . . .  A,  Vi, ref  . .  •  y-^lf){^z2  +  •  •  •  ?/iv)  +  6i (j/i  •  *  •  2/iv V(^XP^-36) 

The  composite  Lyapunov  function  for  the  overall  system  is  given  by:  Vp-i  =  V^_2  + 

y2 

Z_-r^=rl  Ai^p' 

Differentiating  Vp-y  along  the  trajectories  of  the  (xo  •  •  •  Xt,p)  subsystem  and  following 
steps  similar  to  (1.30)-(1.32),  we  obtain 


N 


P-1 


Vp-l  <  -{\min{Qi)  -{p  +  2)}||ei|p  -  -  P*f  -  2'KixQi  Y  xli  -‘^Y 

i=l  _  k=l  2 

T2xz,p{^z,«t— T  Zi^p^i  T  Xi,p—i  T  ^i,p{yi’)  Ai  •  •  •  A,p7  A?  yiyVef  •  •  •  yi^ref) 

+Xi,P^lpiyi,  Zil  ■  ■  ■  kp^  A,  Pi, ref  ■  ■  ■  Piri/)}  + 

+  Y  +  P{d\ik  +  llli/t)}{||Xil|P*’'  +  Ibi.re/IP^}  “  2^*  Qi  Y^  ^xll^  L37) 

k=zl  k~l 

The  following  decentralized  control  input  Ui  is  now  applied  for  the  zth  subsystem 
r  /  r  \^'^i,pXi,p  “h  A,p+1  Xi,p-l  “h  ^t,p(2/n  Zii  .  ,  .  Zi^py  /?i,  yi^ref  •  •  •  ^i^re/) 


Ui  = 


bi^ni~pbi{yi) 

+Xi,„^lf(v„  2.1  ■  2.,(.. ft.  ».,w  .  ■  ■  !/,>./)}  ■ 


(1.38) 


The  properties  of  the  above  designed  decentralized  control  law  are  stated  in  the 
following  theorem. 


Theorem  1.3.1  Suppose  that  the  system  (1.1)  is  transformable  to  the  observer 
canonical  form  (1.2)  and  satisfies  assumptions  (2.2)  —  (2.4).  The  decentralized  control 
(1.38)  along  with  adaptation  law  (1.12)  for  $i  results  in  global  uniform  boundedness  of 
all  the  signals  of  the  closed-loop  system  comprised  of  the  plant  (1.1)  with  disturbance 
a;(t),  the  observer  (1.5),  the  adaptation  gain  Pi  and  the  control  input  ui  (1  <  i  <  N). 
Furthermore,  the  disturbance  signal  does  not  affect  the  tracking  error  (j/j  —  j/j,re/) 
and  the  error  can  be  made  arbitarily  small  based  on  appropriate  choice  of  the  control 
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gains  TTjj.  For  the  regulation  problem,  global  regulation  of  the  states  of  the  plant  and 
observer  is  achieved. 


Proof  of  Theorem  1.3.1:  With  the  control  (1.38),  we  first  show  the  existence  of 
to  counter  the  effects  of  the  interconnections.  One  such  choice  is  given  by 

>  —  max  {dik  +  Vik  +  p{duk  +  yu*)}- 

Qi  k=\...v 

With  \min{Qi)  =  27rio  +  p  +  2,  Xio  >  0,  we  obtain 


N 


t  =  l 


-2xio||ei||^  -  2x,ipi  ^  Xu  -  2  ^ 


k=\  j=2 

2k] 


d"  ^  2^*^{du  +  Vik  +  p{d\ik  +  ^lt7:)}||2/t,re/ 
fc=r 

<  -pV;_i(e,Xi...Xp,/^)  +  0 
where  Xj  =  [xij  •  •  •  XNj]-,  =  [^j . . .  /?;v]  and 


(1.39) 


P  =  ^mm^  (min{2x,oA  J^^(P,),  2xu  . . .  2xip,  Fifr,}) 


(1.40) 


N 


~~  1  X/  +  p{dlik  +  '^lik)}\\yi,ref\ 


2k 


i~l 


k=l 


(1.41) 


'ij^  is  bounded  since  (3*  and  yref  are  bounded.  Thus,  Vp-i{e^Xi  •  —  Xp',  f^)  decreases 
monotonically  along  the  solutions  of  the  (e,  Xi  •  •  •  Xpi  f^)  dynamics  until  it  reaches  the 
compact  set 

fl/  =  {(f,Xi  •  •  -  Xp  J)  e  xn^  x...xTZ^  xTZ^  :  W_i(e,xi  ...XpJ)<Vf^  (1.42) 


where  Vj  =  p  Thus,  the  solutions  (e,  Xi  •  •  •  Xp?  $)  are  globally  ultimately  bounded 
with  respect  to  the  bound  Vj.  The  boundedness  of  other  signals  is  established  as 
follows.  Boundedness  of  Zii(=  yi)  and  rii(xzi,A)  follows  from  that  of  Xzi?  A  and 
yi^ref-  Since  Xi2  is  bounded,  from  (1.11)  we  have  bounded  ^2-  Iteratively,  from 
boundedness  of  yi,  zn.. .  Zim,  yi^rej  •  •  •  Xt',m+i,  and  from  (1.26),  boundedness 

of  A,m+i  is  established  VI  VI  <  m  <  p  ~  1.  Since  ti  are  bounded,  we 

obtain  bounded  zn  . . .  Zip, 

To  prove  boundedness  of  we  use  the  transformation,  rji  —  Tzi^  with 

iVii-^Vip)  =  (^ei**-^tp)  and 


Vij  —  ^ij  ~  7  ^ip-)  ]  —  p  —  1  ... 

^i,Ki  —  p 


(1.43) 
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to  obtain 


Vtj 

Vip 

Vij 

Vi,Ki 


Vi,j+1  +  ■  •  •  Vn),  j  =  l...p -I 

Vi,p-\-l  “I"  “7  Vi,p  "I”  9^j,p(j/l  •  •  •  Vn)  ^i,Ki—p^iiyi)^i 


“1“  Vip 


^i,Ki—j  —  l  ^Z,Avt— p  — 1 


bL-p 


+  { Vi,p+\{yi  ■■■yN)-  -r-^^iAyi  •  •  •  j//v)  ^ ,  i  =  />  + 1 ...«»- 1 


bi,Ki 


^iO 


■^i,p+i 


-  +  I Vi,K,{yi  ■■■yN)  -  -j^^<Pi,p{y^  ■  ■  ■ 

bi,n,-p  I  J 


The  boundedness  of  {Vi,p+i  ^  follows  from  (1.44)  since  rji^p  ~  Zip  and 

(j/i-**yw)  (and  hence  ^ij{yi  •  ^  -  Vn))  are  bounded,  and  the  zero  dynamics  are  ex¬ 
ponentially  stable  (Assumption  (2.2)).  Thus  Zi  =  T~^r]i  are  bounded  V  1  <  i  <  A^. 
Boundedness  of  Ui  now  follows  from  (1.38)  since  Zi^p^i  is  bounded  and  bi^^--p6i{yi)  is 
bounded  away  from  zero.  Boundedness  of  the  states  x  of  the  original  system  (1.1) 
follows  since  2:  =  ^(x)  is  a  global  diffeomorphism. 

The  size  of  residual  set  Clj  given  by  (1.42),  and  hence  the  tracking  error  can  be 
made  arbitrarily  small  by  increasing  the  magnitude  of  the  control  gains  Hij,  Note 
that,  from  definitions  of  //  and  we  have  ij)  (x  ai^  y  (x  cfi  and  oc  F^.  Choosing 
r,-  oc  as  ai  ^  0,  -0  K^ik  +  pdiik)\\yi,rej\\^^ ^  Also,  from  (1.40),  will 

increase  and  the  value  of  n  will  be  determined  by  the  control  gains  Tr^j,  1  <  i  <  /?. 
Thus,  reducing  ai  and  increasing  TTij  will  reduce  the  the  size  of  the  residual  error 
bound  Vf  ~ 

Regulation:  For  the  regulation  problem,  yi^ref  =  0.  In  this  case,  we  set  ai  in  the 
adaptation  law  (1.12)  to  zero.  Thus,  0  =  0.  Following  the  same  procedure  as  in  the 
tracking  case,  Vp-i  is  given  by 


N 


V;>-i(e,Xi---Xp,/5)  <  Y1 

iz=:l 


-2x^0 lie*  11^  -  27r,i^i  ^  xli 


k=l 


(1.45) 


which  implies  that  the  solutions  e,  Xi  •  •  •  Xpi  P  ^re  bounded  for  all  initial  conditions 
and  for  all  t.  A  straightforward  application  of  LaSalle’s  theorem  yields 


lim  ti{t)  =  0,  lim  xpV)  —  0;  1  <  2  <  p;  1  <  «  <  -/V, 

t~^oo  t^OO 

i.e.,  the  tracking  error  xn  —  yi  ~~  yi^ref  and  the  observer  error  Ci  go  to  zero  asymptot¬ 
ically.  □ 

Remark  3.2:  In  contrast  to  [11],  where  p  parameters  need  to  be  adapted  corre¬ 
sponding  to  Y7j-i  0j(pj{y)^  the  adaptive  scheme  proposed  here  requires  adaptation  of 
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only  one  parameter  This  is  true  even  if  the  interconnections  in  (1.2)  were  given  by 

•  •  •  Vn)  long  as  ipij  satisfies  assumption  (1.3).  The  trade-off,  however, 

is  loss  of  asymptotic  tracking  and  utilization  of  high  control  gains  as  in  [12]  to  counter 
the  uncertainties. 

Remark  3,3:  Even  in  the  absence  of  the  disturbance  the  tracking  error  would  be 
bounded,  but  nonzero,  with  a  smaller  size  of  the  residual  set  0/  in  (E42)  since  vuk  = 
Vik  =  0  in  this  case.  The  tracking  error  due  to  disturbance  is  proportional  to  iOmax] 
however,  the  components  (frequencies)  of  the  disturbance  signal  are  asymptotically 
rejected  in  the  output.  With  appropriate  choice  of  control  gains,  this  error  can  be 
made  arbitrarily  small. 

Remark  3.4:  An  important  issue  in  control  of  large  scale  systems  is  to  ensure  that 
the  existing  decentralized  controllers  maintain  stability  and  robust  performance,  if 
subsystems  are  appended  to  the  original  system,  or  taken  off-line,  e.g.,  during  faults 
in  a  power  system.  The  design  methodology  proposed  here  obviates  the  need  for 
controller  redesign  as  long  as  the  order  of  the  interconnections  due  to  the  appended 
subsystem  is  less  than  or  equal  to  that  of  the  original  system.  This  is  generally  true 
for  most  practical  applications  where  the  interconnected  subsystems  are  dynamically 
similar.  In  any  case,  defining  p  as  the  maximum  possible  order  of  all  current  and 
future  interconnections  will  ensure  that  the  same  decentralized  controller  works  for 
the  modified  system. 

Remark  3.5:  From  the  above  design  procedure,  it  is  clear  that  after  the  0th  step,  the 
virtual  controls  can  be  assigned  independently  to  each  subsystem.  Therefore,  even 
the  case  where  the  relative  degrees  of  the  subsystems  are  different,  can  be  treated 
similarly. 

1.4  Conclusion 

In  this  report,  we  have  extended  the  class  of  large  scale  nonlinear  systems  for  which 
decentralized  output  feedback  controllers  can  be  designed.  The  distinguishing  features 
of  this  report  which  makes  it  different  from  all  previously  reported  decentralized 
output  feedback  schemes  are  the  relaxation  of  bounds  on  the  interconnections  from 
linear  to  higher  order  bounds,  the  mismatch  between  the  interconnections  and  the 
control,  and  the  ability  to  reject  any  bounded,  unmeasurable  disturbance  entering  the 
system.  The  stability  and  robustness  properties  attained  are  global.  Although,  the 
scheme  does  not  guarantee  asymptotic  tracking,  the  bounds  on  the  tracking  error  can 
be  made  arbitrarily  small  through  proper  choice  of  control  gains.  For  the  case,  where 
the  objective  is  regulation,  global  asymptotic  regulation  of  all  the  states  of  the  closed- 
loop  system  is  achieved.  One  potential  application  of  the  proposed  scheme  is  in  robust 
decentralized  design  for  large-scale  power  systems  with  swing  angle  measurements 
only  which  is  discussed  in  [14]. 
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2.  Parametrization  of  Stable  Systems  from  Partial  Impulse  Response 

Sequences 

2.1  Introduction 

This  report  addresses  the  problem  of  identifying  the  class  of  all  stable  system 
transfer  functions  that  interpolate  the  given  partial  impulse  response  sequence.  Al¬ 
though  classical  Fade  approximations  match  the  given  impulse  response  sequence  to 
a  maximum  extent  and  are  optimal  in  that  sense,  the  systems  so  obtained  need  not 
be  stable  and  hence  they  may  not  be  attractive  from  physical  considerations.  In  this 
context,  consider  the  problem  of  identifying  a  linear  discrete  time-invariant,  causal, 
stable  system  with  an  unknown  transfer  function  H{z)  from  partial  information  re¬ 
garding  itself.  Since  the  system  is  causal,  it  has  a  one-sided  power  series  expansion 
given  by 

OO 

H{z)  =  Y^hz’^  (2.1) 

k=0 

and  stability  demands  that 

OO 

^\hk\<oo.  (2.2) 

A:=0 

It  follows  from  (2.1),  (2.2)  and  uniform  convergence  that  the  transfer  function  H{z) 
is  analytic  in  \z\  <  1  and  uniformly  continuous^  in  \z\  <  1  [1,2].  Clearly,  the  se¬ 
quence  {hk}'^-Q  represents  the  impulse  response  of  the  system  and  when  the  available 
information  is  of  the  form  k  —  0  ^  the  system  identification  problem  in  the 
rational  case  becomes  equivalent  to  a  Fade  approximation  problem.  In  that  case,  it 
is  easy  to  show  that  rational  ARMA(p,  g)-type  approximations  that  match  the  given 
data  are  unique  provided  p  +  q  <  n  (Fade  approximation)  [1,2].  These  approxima¬ 
tions,  however,  need  not  be  stable  and  hence  from  physical  considerations  they  may 
not  be  acceptable.  For  example,  consider  the  stable  (minimum  phase)  transfer  func¬ 
tion  H{z)  =  .  The  ARMA(1,1)  Fade  approximation  of  this  function  is  given  by 

(2  —  3z)/(2  +  32),  and  it  represents  an  unstable  system  since  Ai{z)  has  a  zero  in 
\z\  <  1. 

In  the  rational  case  the  identification  problem  is  equivalent  to  finding  the  system 
model  order  (p,  q)  and  the  system  parameters.  Given  the  partial  impulse  response  se¬ 
quence,  the  system  model  can  be  established  from  the  invariance  of  the  rank  property 

^Note  that  the  use  of  the  variable  2:  (rather  than  z~^)  here  translates  all  stability  arguments  into 
the  compact  region  \z\  <  1.  H{z)  is  said  to  be  minimum  phase  if  it  is  analytic  together  with  its 
inverse  in  \z\  <  1,  Since  stable  functions  are  free  of  poles  in  \z\  <  1,  in  the  rational  case  they  are 
analytic  in  |2r|  <  1. 
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associated  with  certain  Hankel  matrices  generated  from  this  data.  Thus,  in  particular, 
with  hki  k  >  0,  denoting  its  impulse  response  sequence  as  in  (2.1),  let 

hi  h2  ■■■  hk 
^2  hz  ■  •  •  /ife+l 

hk  hk+1  •  ■  •  h2k--\ 

represent  the  Hankel  matrix  of  size  k  x  k  generated  from  hi,  /i2,  •  •  •  h2k-i-  Then,  for 
a  rational  system  with  degree  p, 

rank  Hk  =  rank  Hp  =  p ,  k  >  p  (2.4) 

and  several  singular  value  decomposition  technique  have  been  proposed  for  model  or¬ 
der  selection  based  on  the  above  rank  condition  [3,4].  Equation  (2.4)  shows  the  linear 
dependence  of  /ip+i,  hp^2^  •  •  •  on  their  p  previous  terms,  and  it  represents  the  finite 
degree  nature  of  a  rational  system.  Although  these  techniques  have  the  advantage 
that  they  can  make  use  of  all  available  impulse  response  data,  they  need  not  lead  to 
stable  systems.  Moreover,  the  above  rank  condition  is  not  valid  in  the  case  of  systems 
that  are  not  rational,  since  they  do  not  represent  finite  degree  systems.  The  problem 
in  that  case  is  to  obtain  equivalent  finite  degree  stable  rational  approximations  that 
capture  all  the  key  features  of  the  original  nonrational  system  in  an  optimal  manner 
by  making  use  of  the  given  data.  Such  a  rational  approximation  should  interpolate 
the  given  information,  and  preferably  be  of  minimum  possible  degree. 

In  this  report,  we  address  this  problem  and  obtain  closed  form  solutions  for  the 
class  of  all  stable  transfer  functions  that  interpolate  the  given  partial  impulse  response 
sequence.  Specifically,  by  making  contact  with  the  Schur  problem  [5]  in,  it  is  shown 
in  section  II  that  the  theory  of  bounded  functions  (Schur  functions)  can  be  utilized 
to  obtain  all  stable  solutions  to  this  problem.  In  this  context,  a  new  model  order 
selection  procedure  is  proposed  here  that  utilizes  the  finite  degree  property  of  a  ratio¬ 
nal  system.  Rational  and  stable  approximation  of  nonrational  systems  is  described  in 
section  III,  by  making  use  of  ideas  developed  in  section  II.  Although  various  authors 
have  addressed  related  problems  in  the  past  utilizing  this  approach  [6]-[14],  some 
interesting  new  observations  will  show  that  rational  system  identification  as  well  as 
stable  rational  approximation  of  nonrational  functions  can  be  realized  from  the  same 
formulation  of  the  Schur  extension  problem. 

2.2  Parametrization  of  Stable  Systems 
2,2.1  The  Schur  Parametrization 

To  start  with,  a  function  d{z)  is  said  to  be  bounded  (Schur  function),  if  (i)  d{z) 
is  analytic  in  \z\  <  1  and  (ii)  \d{z)\  <  1 ,  in  l^j  <  1 .  Because  of  the  analyticity  in 
\z\  <  1,  every  bounded  function  possesses  a  power  series  representation  of  the  form 
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d{z)  =  Yl'k=o  1  kl  <  1 5  that  is  valid  in  \z\  <  1.  If  d{z)  is  rational,  then  \d{z)\  <  1 
in  \z\  <  1  also  implies  d{z)  is  free  of  poles  in  \z\  =  1  and  hence  d{z)  is  analytic  in 
|z|  <  1 .  As  a  result  d{z)  represents  a  stable  system. 

^From  Schur’s  Theorem  [5],  d{z)  defined  above  represents  a  bounded  function  iff 


I-DiDfc>0,  A;  =  0^00,  (2.5) 

where 


do 

0 

0  • 

••  0  ^ 

di 

do 

0  • 

••  0 

d2 

di 

do 

••  0 

(2.6) 

y,  dk 

dfc-i 

dk-2  • 

•  •  do  ) 

represents  the  lower  (or  upper)  triangular  Toeplitz  matrix  generated  from  d,- ,  i  = 

O^k. 

Given  a  partial  set  of  coefficients  ,  A:  =  0  — >  n ,  that  satisfy  I  —  D„D*  >  0, 
“the  problem  of  coefficients”  is  to  obtain  all  bounded  functions  d(z)  such  that  the 
power  series  expansion  of  d{z)  matches  the  given  coefficients,  i.e.,  d{z)  —  dkZ^  + 
0(2"+^) . 

An  algorithm  introduced  by  Schur  in  this  context  answers  this  problem  and  it 
forms  the  basis  for  our  approach  to  the  present  parametrization  problem.  As  Schur 
has  first  observed,  if  d{z)  represents  a  bounded  function,  then,  so  does  the  function[5] 


di(^)  = 


1 

2: 


d{z)  —  d(0) 

1  —  d*(0)d(z)  ’ 


d(0)  =  do 


(2.7) 


In  the  rational  case,  since  z  =  0  is  not  a  pole  of  di(z),  from  (2.7)  we  obtain  that 
the  degree^  of  the  new  bounded  function  di(z)  never  exceeds  that  of  d(z),  i.e.. 


d(di(z))<d(d(^)),  (2.8) 

with  inequality  iff  the  Ijz  factor  in  (2.7)  cancels  a  pole  of  (d(2:)  — do)/(l  —  dod(2r))  [17]. 
Since  this  cancellation  can  occur  only  at  z  =  oo,  from  (2.7)-(2.8),  degree  reduction 
happens  iff  the  denominator  term  1  —  d(0)*d(2)  satisfies  1  —  d*(0)d(2;)|^_^  =  0 ,  or 

S(di(z))  <  S(d(z))  4^  d(z)d*(2)|^^o  =  1 ,  (2.9) 

where 

d,(z)  =  d*(l/2*)  (2.10) 


^The  degree  S(II(z))  of  a  rational  function  B'(z)  equals  the  totality  of  its  poles  (or  zeros),  with 
multiplicities  counted,  including  those  at  infinity. 
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is  defined  to  che  paraconjugate  form  of  d{z).  Clearly,  the  paraconjugate  form 
represents  the  ordinary  complex  conjugate  operation  on  the  unit  circle. 

The  bilinear  transformation  in  (2.7)  maps  the  inside  of  the  unit  circle  onto  itself. 
Thus,  in  general,  with  dk{z)  representing  an  arbitrary  bounded  function  and  with 


Sk  =  4(0) , 

(2.7)  translates  into 

j  (  \  ^  dk(^Z^  Sk  ;  \  n 

z  1-Sf.dk{z) 


(2.11) 


(2.12) 


with  the  understanding  that  do{z)  =  d[z).  The  above  Schur  algorithm  can  be  recur¬ 
sively  updated  and  after  n  such  steps,  we  get 


\  K{z)  zdn{z)dnJ,\{z) 

(l\^Z  j  ^  ^ 

a„(z)  -f  zbn{z)dr,^-i{z) 

where  an(z)  and  6„(z)  are  in  general  two  polynomials  of  degree  n,  and 

dn(z)  =  z'^an*{z)  =  2:"a*(l/2:*)  , 


(2.13) 


(2.14) 


k{z)  =  z^bUz)  =  z^hl{llz*) 


(2.15) 


represent  polynomials  reciprocal  to  an[z)  and  bn{z)  respectively  that  satisfy  the  re¬ 
cursion 

—  ^71— 1('^)  d~  Z ,  77/  ^  1  ^2.16^ 

and 

bn{z:)  =  hn-l{z)  +  ZSnan-\{z)  ,  72  >  1.  (2.17) 

an{z)  and  bn{z)  are  defined  to  be  the  Schur  polynomials  of  the  first  and  second  kind 
respectively.  Notice  that,  if  d{z)  is  rational  to  start  with,  application  of  the  above 
procedure  will  result  in  a  rational  bounded  function  dn-\-i{z)  for  every  rz,  and,  further 
from  (2.8),  in  that  case  6{dn-\-i{z))  <  6{d{z))  ^  n  >0.  ^From  (2.13),  the  iterations 
in  (2.16)“(2.17)  start  with 

ao(^)  =  1 ,  bo{z)  =  So  =  do.  (2.18) 

Using  (2.16)-(2.17),  it  is  easy  to  show  that  an{z)  ^  n  =  1  — >  cxd,  represent  strict 
Hurwitz  polynomials^.  In  particular  with  dn-\-i{z)  =  0 ,  we  have  d{z)  =  bn{z)fan{z) 

Hurwitz  polynomial  is  free  of  zeros  in  \z\  <  1,  and  a  strict  Hurwitz  polynomial  is  free  of  zeros 
in  hi  <  1. 
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itself  is  bounded  and,  moreover,  a  direct  expansion  gives  d{z)  — 

i.e.,  the  power  series  expansions  of  the  bounded  functions  d{z)  and  hn{z)  j an{z)  agree 

upto  the  first  n  +  1  terms  and  hence  the  above  terms  must  be  independent  of  {z) , 

Thus,  for  every  arbitrary  bounded  function  <i„+i(2),  we  must  have  the  interpolation 

property 

=  K{Z)  +  Zdr.{z)dn+x{z)  ^  ^ 

+  to 

and  the  d^s,  k  =  0  ^  n,  can  be  determined  from  the  Schur  polynomials  0^(2), 
bk{z)  in  (2.16)-(2.18).  Conversely,  (2.19)  is  completely  specified  by  the  first  (n  +  1) 
coefficients  {dk}k=oi  or,  from  the  Schur  polynomials  a„(2)  and  bn{z).  To  complete 
the  recursions  in  (2.16)-(2.18),  only  the  coefficients  ^  =  0  — >  n,  are  required  and 
they  can  be  obtained  recursively  from  the  given  data  dk,  k  —  0  n  as 

So  =  do  (2.20) 


and 


(2. 


where  {  }„  represents  the  coefficient  of  in  {  }.  Using  this,  an{z)  and  bn{z)  can 
be  computed  recursively,  and  the  class  of  all  bounded  functions  that  interporate  the 
given  coefficients  dk,  k  =  0  n  is  given  by  d{z)  in  (2.19). 

In  general,  the  given  impulse  response  data  hk ,  k  =  0  n ,  do  not  form  part 
of  a  bounded  function,  and  to  make  use  of  the  above  formulation  in  section  II,  it 
is  necessary  to  ‘prepare’  this  data  so  that  it  confirms  with  a  bounded  function.  To 
attain  this  goal,  consider  the  matrix 


— 


^4”  ^Un-k 

k=0 


k=-l 


n—1 


1  - 


k=0 


1-  (4-1(2)^424 

t  ^=0  J  n— 1 


,  n  >  1 


H„  = 


/ho  0  •••  0  \ 

hi  ho  •  •  •  0 

\  hji  hn-i  ■  •  •  ho  / 


(2.22) 


and  let  Aj(n)  represent  the  largest  eigenvalue  of  H„H*.  Then,  clearly,  the  sequence 
dk  =  —  ,  Kn  >  Ai(n) ,  /c  =  0  — >  n  ,  (2.23) 
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satisfies  (2.5)  with  inequality,  and  hence  qualifies  as  the  first  n  +  1  coefficients  of  a 
bounded  function.  Recursive  determination  of  the  coefficients  Sk ,  k  =  0  n  from 
(2.20)-(2.21),  together  with  ak{z)  and  bk{z),  k  =  0  n  using  (2.16)-(2.18),  gives 


H{z)  =  Kn 


K{z)  +  Zan{z)dn+l{z) 
an{z)  +  zbn{z)dn^i{z) 


Y,hkz’^^0{z^^^) 

k=0 


(2.24) 


to  be  the  class  of  all  transfer  functions  that  are  analytic  in  l^j  <  1,  free  of  poles  in 
\z\  <  1  and  interpolate  the  given  partial  impulse  response  sequence  hk,  A:  =  0  — >  n . 

Equation  (2.24)  can  be  given  two  interesting  interpretations:  First,  if  a  system 
transfer  function  H{z)  is  rational  to  start  with,  then  its  representation  as  in  (2.24) 
after  n  steps  of  the  Schur  algorithm  will  imply  that  dnj^\{z)  must  be  a  rational  func¬ 
tion.  Similiarly  if  H{z)  is  nonrational  to  start  with,  then  dn+i{z)  must  be  nonrational 
in  (2.24). 

The  alternate  interpretation  shows  that  given  /io,  hi,  ■  ■  ■  /i„,  equation  (2.24)  rep¬ 
resents  all  stable  system  transfer  functions  both  rational  and  nonrational  that  inter¬ 
polate  the  given  data,  and  they  can  be  obtained  by  varying  dn+i{z)  over  all  bounded 
functions.  Thus  even  if  the  given  data  corresponds  to  a  nonrational  system,  the 
freedom  present  in  the  choice  of  dn^i{z)  in  (2.24)  can  be  utilized  for  rational  approx¬ 
imation  of  H{z)  by  appropriate  choice  of  rational  bounded  functions  dn^i{z). 

The  above  discussion  shows  that  d„_|_i(0)  can  be  utilized  for  rational  system  iden¬ 
tification  as  well  as  rational  approximation  of  nonrational  systems.  In  particular,  if 
dn+i{z)  is  chosen  to  be  a  rational  bounded  function,  then  since  dn+i{z)  and  H{z)  are 
free  of  poles  in  \z\  <  1,  H{z)  in  (2.24)  represents  a  stable  regular  rational  transfer 
function  (analytic  in  \z\  <  1)  that  matches  the  given  coefficients.  As  a  result,  the  class 
of  all  stable  rational  functions  that  interpolate  the  given  impulse  response  sequence 
is  obtained  from  (2.24)  by  varying  dn+i{z)  over  all  rational  bounded  functions. 


2.2.2  The  Rational  Case 

If  d{z)  in  (2.13)  is  rational  to  start  with,  repeated  application  of  the  Schur  pro¬ 
cedure  will  result  in  rational  bounded  functions  d„^i[z)  that  satisfy  the  degree  con¬ 
straint  in  (2.8)-(2.9).  As  a  result,  from  (2.19)  and  (2.24),  it  follows  that  every  stable 
rational  function  H{z)  can  be  represented  as  in  (2.24)  where  d„+i(2r)  is  a  unique 
rational  bounded  function  that  satisfies 

8{d^^i{z))<8{H{z)),  (2.25) 

and  degree  reduction  in  (2.25)  happens  according  to  (2.9).  Thus  if  H{z)  represents  a 
stable  ARMA(p,  q)  system  with  p  >  ^  in  (2.24),  then  8{H{z))  =  p  and  since  =  0 , 
it  follows  from  that 

8{dn+i{z))=p  (2.26) 
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and  further  using  the  degree  arguments  as  in  (2.8),  we  obtain 


/o  +  /i-2^  +  •  •  •  +  * 

1  +  giz  +  g2Z  +  •  •  •  +  gpZ^  ’ 


n  >  1 . 


(2.27) 


(Since  g{z)  is  strict  Hurwitz,  go  ^  0  and  it  is  normalized  here  to  unity.)  Substituting 
(2.27)  into  (2.24)  we  get 


H{z)  = 


bniz)g{z)  +  zf{z)an{z) 


an{z)g{z)  +  zf(z)bn(z)  k=0 


=  '£hkZ^  +  0{z-^^). 


(2.28) 


Since  every  rational  system  after  repeated  application  of  the  Schur  procedure  has  the 
above  representation  for  any  n,  where  f{z)fg{z)  is  a  unique  rational  bounded  function 
as  in  (2.27),  we  can  make  use  of  the  degree  constraint  of  H(z)  in  (2.28)  to  obtain 
this  unknown  bounded  function.  Towards  this,  notice  that  the  formal  degree  of  both 
the  numerator  and  denominator  of  (2.28)  is  n  +  p,  and  to  respect  the  ARMA(p,  9) 
nature  of  H(z),  we  first  equate  the  denominator  coefficients  of  z^'^'^,  •  •  • , 
to  zero.  However,  as  shown  in  [18],  equating  the  coefficients  of  •  •  • ,  2^+" 

in  the  denominator  to  zero  implies  that  the  respective  coefficients  in  the  numerator 
are  also  zeros.  As  a  result,  we  obtain  n  equations  from  the  denominator  coefficients  of 
^p+i^  _-,p+2^  . . .  ^  ^p+n  ^  p  —  q  equations  from  the  remaining  numerator  coefficients 
of  •  •  • ,  z^  .  Thus  we  have  n  p  —  q  equations  and  2p  unknowns  gk  ,  k  =  I  p 
and  fk,  k  =  0  p  —  1 .  Clearly  the  minimum  value  of  n  is  given  hy  n  =  p  q  and 
in  that  case  the  resulting  2p  equations  in  2p  unknowns  can  be  represented  in  matrix 
form  as 


Ax  =  b 


(2.29) 


where 


^p-hq 

0  •• 

0 

0 

bo 

0 

••  0 

0 

^p+q—l 

O-p+q 

0 

0 

bl 

b*o  ■ 

••  0 

0 

aq^2 

O-g+3 

Gp+g 

0 

K-2 

bl-3 

bo 

0 

aq+2  ■  ■ 

^p+g— 1 

Gp^q 

K-r 

K-2  • 

■■  bl 

b*o 

Gi 

02  •  • 

Gp-i 

Gp 

K+1-i 

h* 

^p+g-2 

K 

bo 

6i  •• 

bp-2 

bp—i 

K+i 

<+l 

0 

bo 

bp— 3 

CA 

1  . . 

0 

«g+2 

0 

0  •• 

•  b, 

bq+1 

0 

0 

«p-i 

0 

0  •• 

bq  —  l 

bq 

0 

0 

9i  1  /p-i 

fp-2 

•••  /l 

gp  9p-i 

... 

/o  ] 

(2.30) 


(2.31) 


and 


0  0 


0  ^p-\-q 


^p— 


1  .  .  . 


(2.32) 


Here  ak^bk^k=^0^p  +  q  represent  the  coefficients  of  the  p  +  q  degree  Schur 
polynomials  ap^q{z)  and  bp^q{z)  respectively,  i.e., 


Gp^q^z)  —  Go  +  GiZ  +  a2Z  Gp^gZ^'^^  , 


(2.33) 


and 


bp^q[z)  —  bo  biZ  b2Z  bp^qZ^^^  . 


(2.34) 


Note  that  p  q  represents  the  minimum  value  for  the  available  data  points  n,  and 
if  a  larger  number  of  such  data  is  available,  then  these  equations  can  be  modified  to 
accomodate  that  situation  leading  to  an  overdetermined  system  of  equations  in  (2.29). 
At  the  correct  stage,  (2,29)  is  guaranteed  to  have  a  unique  solution  that  results  in  a 
bounded  function  for  f{z)/g{z)^  and  the  unknown  system  parameters  of  H{z)  can  be 
expressed  in  terms  of  the  gk^s  and  fk'^s  so  obtained.  In  fact,  from  (2.28),  with 


H{z) 


m 

P(z) 


Qo  +  QlZ  +•••-!-  QqZ^ 

Pq  +  P^Z  +  *  •  *  +  PpZ^ 


k=0 


(2.35) 


we  get 
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(2.36) 


k  k—1 

Pk  =  Y^  aigk-i  +  K-ifk-i-i ,  k  =  0^p 


i=0 

t=0 

and 

Qk  — 

/  k  k-l  \ 

(  X]  ^i9k~i  +  1 

\i=0  i-0  / 

T 

O 

11 

(2.37) 

Clearly,  stability  of  H{z)  and  the  interpolation  property  follows  from  the  bounded 
character  of  dn+i{z)^  and  since  P{z)  and  Q{z)  are  computed  without  involving  any 
spectral  factorization,  the  nonminimum  phase  characteristics  of  H{z)  if  any  is  also 
preserved  here. 

Finally,  to  show  the  uniqueness  of  the  ARMA(p,  5)  form  in  (2.35),  the  degree 
reduction  possibility  through  common  factor  cancellation  in  its  numerator  and  de¬ 
nominator  must  be  ruled  out.  To  see  this,  let  dn+-i.{z)  =  f{z)fg{z)  represent  any 
rational  bounded  function  in  (2.28).  Then  from  (2.19)  and  (2.28),  since  the  denomi¬ 
nator  P{z)  =  an{z]  -h  zdn+i{z)bn{z)  is  strict  Hurwitz,  any  such  common  zero  between 
P{z)  and  Q{z)  =  bn{z)  zdn+i{z)dn{z)  must  be  outside  the  unit  circle.  Let  zq  repre¬ 
sent  such  a  common  zero.  Then  |.2ro|  >  1  >  and  from  P{zo)  =  Q(zo)  =  0,  it  is  easy  to 
show  that 

|d„+i(^o)|>l.  (2.38) 

On  the  other  hand,  referring  back  to  P(zo)  =  0 ,  we  have 

zod„+i(zo)^-^4H,  (2.39) 

and  since  an(z)  is  free  of  zeros  in  |^|  >  1,  the  function  bn(z)/an(z)  is  analytic  in 
1^:1  >  1  and  it  is  bounded  by  unity  on  the  unit  circle.  Thus  by  maximum  modulus 
[15],  bn(z)/an(z)  is  bounded  by  unity  everywhere  in  \z\  >  1,  and  in  particular  from 
(2.39),  at  2:  ==  we  have 

l^odn+i(2:o)|  =  -T-y  \  <  1, 

or  |d„+i(zo)|  <  1,  which  contradicts  (2.38),  implying  that  the  numerator  and 
denominator  factors  in  (2.35)  are  free  of  any  such  common  factors  outside  the  unit 
circle,  and  hence  in  the  entire  2-plane.  This  proves  the  uniqueness  of  the  degree 
reduction  condition  through  the  procedure  described  in  (2.29)-(2.32). 

Equations  (2.27)-(2.37)  can  be  implemented  provided  p  and  q  are  known.  Usu¬ 
ally  the  model  order  (p,  q)  is  unknown,  and  that  will  have  to  be  evaluated  from  the 
given  data.  As  we  show  below,  the  invariant  characteristics  of  the  rational  bounded 
function  dn+i{z)  in  (2.26),  together  with  the  Schur  update  rule  in  (2.12)  can  be  used 
to  determine  the  model  order. 
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2,2.3  Model  Order  Selection 


Having  determined  ^^+,+1(2)  =  f{z)lg{z)  as  in  (2.29)-(2.32),  the  bounded  func¬ 
tion  dp^q^2{z)  at  the  next  stage  [n  =  •p  +  q-\-\)  can  be  evaluated  in  a  similar  manner 
from  the  Schur  polynomials  ap+q+i{z)  and  bp+q+i{z).  In  fact,  letting 


from  (2.28),  we  also  have 

^  ^  .  6p+g+i(z)e(^)  -I-  zdp+q+i{z)c{z) 

ap+q+iiz)e{z)  -h  zbp+q+i{z)c{z) 


(2.40) 


(2.41) 


and  as  before  e(z)  and  c{z)  can  be  evaluated  by  equating  the  coefficients  of 

•  •  •,  in  the  denominator  and  •  •  • ,  z^  in  the  numerator  to 

zero'*.  Once  again,  these  equations  possess  a  unique  solution  at  the  correct  stage  for 
the  unknowns  ek.,k  =  l^p,Ck,k  =  0^p— 1,  (with  cq  =  1),  and  they  result  in  a 
bounded  function  in  (2.40).  Notice  that  both  these  bounded  functions  dp+q+i(z)  and 
dp^q+2{z)  are  of  degree  p,  have  the  same  form  as  in  (2.27),  and  are  related  through  the 
Schur  rule  as  in  (2.12).  dp+q+2{z)  makes  use  of  additional  information  hp+g+i  about 
the  system,  through  the  new  Schur  polynomials.  Substituting  these  two  bounded 
functions  into  (2.12)  and  simplifying,  we  obtain 


f{z)  _  foe{z)  -I-  zc{z) 
g{z)  e(z)  +  zf^c(z)  ' 


(2.42) 


Equation  (2.42)  relates  the  coefficients  of  the  bounded  functions  at  two  consecutive 
stages,  and  equating  the  ratios  of  like  powers  on  both  sides  of  (2.42)  and  rearranging, 
we  obtain  the  conditions 


^k(p,g)  =  0,  k  no¬ 
where 

^o(Pi  g)  =  fo^p  +  Cp_i 

and 

A  fo^k  +  Cfc-l 


P-  1, 


^k{p,q) 


+  fo^k- 


h 

9k 


k  =  1 


1. 


(2.43) 


(2.44) 


'^Although  this  results  in  {2p+  1)  equations  in  2p  unknowns,  since  the  coefficient  of  is  the 

same  as  (2.43),  the  remaining  2p  equations  are  implemented  in  our  computations. 
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These  conditions  are  a  direct  consequence  of  (2.26),  and  reflect  the  AR]ViA(p,  q)  nature 
of  the  problem.  Since  the  first  stage  where  (2.43)  and  (2.44)  are  satisfied  occurs  at 
the  correct  stage,  by  updating  p  and  q  sequentially  starting  with  p  >  1  and  q  <  p, 
the  true  model  order  can  be  found  as  the  smallest  integers  p  and  q  that  satisfy 


^o{p,q)  =  0, 


or,  more  genreally 


e(p,  q) 


m 


'I 


Y.  =  0- 


(2.45) 


The  key  feature  of  a  rational  system  -  its  degree  -  is  exploited  here  in  determining 
the  true  model  order  and  system  parameters. 

The  nonminimum  phase  stable  rational  transfer  function  examples  in  Figs.  1-2 
highlight  all  important  aspects  of  the  algorithm  described  in  this  report.  The  unknown 
system  is  assumed  to  be  ARMA(n,  m)  with  n  >  m,  and  initialization  begins  with  n  — 
1,  m  =  0.  Preparation  of  the  given  impulse  response  sequence  hk-,  A:  =  0— >-n  +  m  +  l, 
is  first  carried  out  to  generate  dk  =  hkj  Kn+m+ii  /i;  =  0-^n  +  m  +  lasin  (2.22)-(2.23). 
Computation  of  the  Schur  polynomials  a„+TO+i(2),  bn+m+i{z)  using  (2.16)-(2.17)  and 
(2.20)-(2.21),  followed  by  those  of  the  bounded  functions  dn+m+i{z)  =  f{z)lg{z)  and 
dn^rn+2{z)  =  c(z) ( e(z)  then  allow  eo(n,m)  and  e{n,m)  to  be  evaluated  using  (2.43), 
(2.45),  provided  both  dn+m+i{z)  and  dn+rn+2{z)  exist  as  bounded  functions.  The  heavy 
dots  on  all  curves  in  Figs.  l(c)-2(c)  indicate  the  presence  of  such  a  stage  (n,m),  and 
if  such  is  not  the  case  that  particular  stage  is  skipped  and  the  indices  n  and  m  are 
updated.  Notice  that  eo{n,  m)  and  e(n,  m)  are  guaranteed  to  exist  at  the  correct  stage 
n  =  p  and  m  =  q,  and  since  the  first  place  where  eo(n,  m)  and  e(n,  m)  equal  zero  also 
occurs  at  the  correct  stage,  sequential  updating  of  n  and  m  continues  until  substantial 
relative  minima  in  the  values  of  eo{n,m)  and  e{n,m)  are  observed  to  occur  for  the 
first  time.  The  corresponding  pair  (n,m)  is  then  identified  as  the  model  order  {p,q) 
and  the  system  parameters  are  computed  from  (2.35)-(2.37).  Finally,  to  facilitate 
comparison,  the  exact  magnitude  \H{e^^)\  and  its  reconstructed  counterpart  (dotted) 
are  plotted  in  Figs.  l(a)-2(a).  Similarly,  the  exact  phase  (j){9)  and  the  reconstructed 
phase  (dotted)  are  plotted  in  Figs.  l(b)-2(b). 

Although  the  theoretical  development  in  section  III  assumes  p  >  as  the  MA(5) 
example  in  Fig.  2  shows,  every  case  where  q  >  p  can  be  detected  as  an  ARMA(g,9) 
system.  This  means  of  course  that  some  of  the  reconstructed  coefficients  in  the 
denominator  are  filled  in  automatically  as  zeros,  to  raise  the  denominator  degree  to 
<1- 
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2.3  Stable  Rational  Approximation  of  Nonrational  Systems 

A  nonrational  system  has  a  transfer  function  that,  unlike  the  rational  systems, 
cannot  be  expressed  as  the  ratio  of  two  polynomials  of  finite  degree.  If  such  a  system 
is  stable,  then  it  admits  a  power  series  expansion  in  \z\  <  1,  and  the  problem  is  to 
represent  this  by  a  rational  system  in  some  optimal  manner.  As  remarked  in  the 
introduction,  although  Fade  approximations  can  achieve  this  goal,  such  approxima¬ 
tions  need  not  guarantee  stability.  Moverover,  the  Hankel  matrices  generated  from 
the  impulse  response  data  has  no  particular  rank  invariant  structure  in  this  case. 
In  this  context,  once  again  we  can  make  use  of  (2.24)  to  obtain  all  stable  rational 
solutions  to  this  problem. 

To  start  with  notice  that  the  Schur  extraction  principle  in  (2.12)-(2.13)  is  per¬ 
fectly  general  and  hence  if  H{z)  is  nonrational,  then  the  representation  in  (2.24) 
is  still  valid,  where  d„+i(z)  in  that  case  represents  a  nonrational  bounded  function. 
Interestingly,  as  remarked  there,  if  dn+i{z)  is  replaced  by  any  rational  bounded  func¬ 
tion,  we  obtain  a  stable  rational  transfer  function  Hr{z)  that  interpolates  the  given 
impulse  response  sequence  This  key  observation  can  be  used  to  determine 

interpolating  rational  systems  with  minimum  degree. 

Since  Fade  approximations  preserve  the  optimal  degree  character,  if  such  ap¬ 
proximations  are  also  stable,  then  they  must  follow  from  (2.24)  for  a  specific  ratio¬ 
nal  bounded  function  ^^+1(2).  To  determine  such  bounded  functions,  let  dn+i{z)  = 
f{z)fg{z)  represent  a  degree  m  bounded  function  that  when  substituted  into  (2.24) 
generates  an  ARMA(p,  g)  Fade-approximation  Hr{z).  Thus 


Hr{z) 


Q{z)  _  _  K{z)g{z)  +  zdn{z)f{z) 

P{z)  ”  ar,{z)g{z)  +  zbn{z)f{z) 


J2hz'^  +  0{z^+^). 

kz=:0 


(2.46) 


For  (2.46)  to  represent  the  Fade  approximation,  we  must  have  p  +  9  <  n,  and  once 
again  to  respect  the  ARMA(p,  g)  nature  of  Hr{z)^  the  polynomial  f{z)  must  be  of 
degree  m  —  1  and  hence  g{z)  must  be  of  degree  m.  Thus  the  formal  degrees  of  P{z) 
and  Q{z)  in  (2.46)  are  n  -f  m,  provided  8[an{z))  =  n,  and  hence  the  coefficients 
of  •  •  • ,  2”+’"  in  the  denominator,  and  the  coefficients  of  •  •  • ,  2”+™  in  the 

numerator  must  be  zeros.  It  can  be  shown  that  the  coefficients  of  2:”+'" 

in  the  numerator  and  denominator  generate  the  same  equations  and  hence  this  gives 
n  + {m  —  p)  +  [m  —  q)  =  n  +  2m  —  {p  +  q)  equations  in  2m  unknowns.  Since  n  >  p  +  q, 
there  are  at  least  2m  equations  and  they  can  be  used  to  solve  for  the  unknowns.  ^From 
the  above  degree  argument  m  —  p  >  0  and  m  —  5  >  0,  or  m  >  max(p,  q),  and  hence  for 
a  given  p,  q  (with  p'>  q),  the  least  complex  bounded  function  dn+i{z)  is  also  of  degree 
p.  In  that  case,  the  desired  bounded  function  dn+i{z)  has  exactly  the  same  form  as  in 
(2.27),  and  the  system  of  equations  so  obtained  has  the  functional  representation  in 
(2.29)-(2.32).  However,  unlike  the  rational  case,  the  system  of  equations  so  obtained 
need  not  yield  a  solution  for  g{z)  and  f{z),  and  even  if  a  solution  exists  there,  g{z)  so 
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obtained  need  not  be  strict  Hurwitz,  and  further  f[z)fg{z)  need  not  turn  out  to  be 
a  bounded  function.  However,  for  some  p,  q,  if  f{z)lg{z)  turns  out  to  be  a  bounded 
function,  then  Hr{z)  in  (2.46)  represents  a  stable  ARMA(p,  5)  Fade  approximation  to 
the  given  nonrational  function.  Thus  every  stable  Fade  approximation  to  the  given 
data  has  the  representation 


R  (A  =  ^  bp+g{^)g{z)  +  zap+^{z)f{z) 

P(z)  ap+,{z)g{z)  +  zbp+g(z)f{z) 


^  hkz'^  +  Oiz^+^) , 


k~0 


(2.47) 


where  n  >  p  +  q,  and  f{z)(g{z)  represents  a  bounded  function  given  by  (2.29)-(2.32). 
We  summarize  the  above  observations  as  follows: 


The  necessary  and  sufficient  condition  for  the  existence  of  a  stable  ARMA(p,  9)  Fade 
approximation  to  the  impulse  response  sequence  {hkYkto  system  of  linear 

equations  in  (2.29)-(2.32)  generated  from  the  associated  Schur  polynomial  coefficients 
yield  a  bounded  solution  of  degree  q  for  f(z)fg{z)  in  (2.27).  In  that  case,  (2.35)-(2.37) 
and  (2.47)  represent  the  desired  stable  transfer  function. 

Interestingly,  the  above  remarks  raise  the  following  question:  Given  an  H{z)  that 
represents  a  stable  (nonrational)  system  transfer  function,  does  there  always  exist  a 
stable  ARMA(p,  q)  Fade  approximation  for  some  p  and  Clearly,  if  such  a  solution 
exists,  then  that  must  follow  from  (2.47)  with  n  >  p  +  q  for  a  rational  bounded 
function  f{z)fg{z)  of  degree  p,  that  is  obtained  by  solving  the  system  of  equations  in 
(2.29)-(2.32). 

Although  it  is  easy  to  show  that  the  strict  Hurwitz  character  of  g{z)  must  occur 
infinitely  often,  indicating  the  possibility  of  the  desired  bounded  function  in  (2.46)- 
(2.47),  unfortunately,  nothing  more  specific  can  be  said  about  the  bounded  character 
of  d„.|_i(.2:)  J {^z'^ j g^z'^t  In  fact,  given  k  —  0  ^  n,  there  might  exist  no  stable 
Fade  approximation,  and,  moreover,  the  lowest  nontrivial  stable  rational  approxima¬ 
tion  that  matches  the  given  data  might  be  of  ARMA(n,n)  form  which  corresponds 
to  dn+i  =  0  in  (2.24).  However,  so  long  as  f{z)/g{z)  is  chosen  to  be  a  rational 
bounded  function  in  (2.46),  it  represents  the  class  of  all  stable  rational  transfer  func¬ 
tions  that  interpolate  ho,  hi,  ■  ■  ■  hn-  Thus  even  if  stable  Fade  approximations  and 
absent  in  a  particular  situation,  by  relaxing  the  Fade  constraint,  other  stable  rational 
approximations  can  be  obtained. 

Figs.  3-4  as  well  as  our  extensive  computations  involving  nonrational  systems 
containing  logarithmic  and  essential  singularities  seem  to  indicate  that  stable  Fade 
approximations  always  exist.  In  the  simulations  pi'esented  here,  the  original  nonra¬ 
tional  function  H{z)  is  used  to  compute  h^,  k  —  ^  ^  p  q,  and  it  is  first  ‘prepared’ 
to  generate  djt.  A:  =  0  — >  p  -f  g,  as  in  (2.22)-(2.23)  and  thereby  the  Schur  poly¬ 
nomials  Gp+,(z)  and  hpj,.g{z).  From  the  above  theorem,  since  a  stable  ARMA(p,  9) 
Fade  approximation  to  this  data  must  follow  from  (2.47)  for  a  bounded  solution  of 
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f{z)lg{z)  given  by  (2.29)-(2.32),  those  equations  are  verified  for  such  a  solution,  and 
the  indices  are  updated  in  a  sequential  manner.  The  heavy  dots  in  Figs  3(c)-4(c) 
indicate  the  presence  of  such  a  stage,  where  (2.29)-(2.32)  yield  a  bounded  solution 
for  =  f{z)fg{z)  that  results  in  a  stable  Fade  approximation.  Further,  if  the 

two  consecutive  rational  functions  dp^q^i{z)  and  dpj^q^2{^)  turn  out  to  be  bounded 
and  are  related  as  in  the  Schur  algorithm,  then  eo{p,q),  ^k{p,q)  will  be  zero.  When 
H{z)  is  known  in  advance,  the  percent  spectral  error 


r}{p,q)  =  sup 
6 


(2.48) 


also  may  be  used  for  model  order  selection.  Notice  that  Fig.  4  represents  a  transcen¬ 
dental  nonminimum  phase  system  (zero  at  the  origin)  with  a  logarithmic  singularity 
at  z  =  1.  Nevertheless,  as  seen  from  Figs.  4(a)-(c),  the  ARMA(14, 4)  Fade  ap¬ 
proximation  is  stable  and  preserves  the  nonminimum  phase  character  of  the  original 
system.  The  abundance  of  stable  Fade  approximations  are  evident  in  Figs  4(c)-5(c). 
Nevertheless  a  rigorous  proof  is  still  lacking  in  the  general  case  regarding  the  bounded 
character  of  f{z)/g{z)  for  some  p,  q,  and  the  issue  remains  unresolved. 

2.4  Conclusions 

This  report  investigates  the  problem  of  obtaining  all  stable  rational  solutions  that 
interpolate  the  given  partial  impulse  response  sequence  by  making  use  of  the  well 
known  theory  of  bounded  (Schur)  functions.  In  this  context,  a  new  test  criterion  is 
developed  to  determine  the  model  order  of  rational  systems,  and  thereby  determine 
their  system  parameters  from  the  given  impulse  response  sequence.  The  theory  de¬ 
veloped  is  further  utilized  to  obtain  the  necessary  and  suflficient  conditions  for  stable 
Fade  approximations  of  nonrational  systems.  A  practical  algorithm  is  developed  that 
translates  the  stability  condition  into  the  bounded  character  of  a  rational  function 
generated  from  a  set  of  linear  equations  obtained  from  the  Schur  polynomial  coeffi¬ 
cients  associated  with  the  given  impulse  response  sequence.  Interestingly,  since  the 
present  technique  does  not  make  use  of  any  factorization  procedure,  the  nonminimum 
phase  characteristics  of  the  original  system  are  preserved  here. 


2.5  References 

[1]  G.  A.  Baker,  Jr.,  Essentials  of  Fade  Approximants,  New  York:  Academic  Fress, 
1975. 

[2]  W.  B.  Jones  and  W.  J.  Thron,  Continued  Fractions:  Analytic  Theory  and  Ap¬ 
plications^  Cambridge  University  Fress,  1984. 


31 


[3]  Y.  T.  Chan  and  J.  C.  Wood,  “A  new  order  determination  technique  for  ARMA 
processes,”  IEEE  Trans.  Acoust.,  Speech,  Signal  Processing,  vol.  ASSP-32,  no.  3, 
pp.  517-521,  June  1984. 

[4]  J.  J.  Fuchs,  “ARMA  order  estimation  via  matrix  perturbation  theory,”  IEEE 
Trans.  Autom.  Control,  vol.  AC-32,  no.  4,  pp.  358-361,  April  1987. 

[5]  I.  Schur,  “On  power  series  which  are  bounded  in  the  interior  of  the  unit  circle,” 
Parts  I  and  II,  in  I.  Schur  Methods  in  operator  theory  and  signal  processing  ed. 
by  I  Gohberg,  pp. 31-88,  Birkhauser,  Boston,  1986. 

[6]  H.  Dym,  ^^J-Contractive  Matrix  Functions,  Reproducing  Kernel  Hilbert  Spaces 
and  Interpolation,  No.71  in  CBMS  regional  conference  series.  Providence:  Amer¬ 
ican  Math.  Soc.,  1989. 

[7]  P.  Dewilde  and  H.  Dym,  “Schur  Recursions,  Error  Formulars,  and  Convergence 
of  Rational  Estimators  for  Stationary  Stochastic  Sequences,”  IEEE  Trans.  In¬ 
format.  Th.,  vol. 27,  pp446-461,  July  1981. 

[8]  Ph.  Delsarte,  Y.  Genin,  and  Y.  Kamp,  “On  the  role  of  the  Nevanlinna-Pick 
Problem  in  Circuit  Theory  and  Design,”  Circuit  Th.  AppL,  vol. 9,  ppl77-187, 
1981. 

[9]  D.  Alpay,  P.  Dewilde,  and  H.  Dym,  “On  the  Existance  and  Construction  of 
Solutions  to  the  Partial  Lossless  Inverse  Scattering  Problem  with  Applications 
to  Estimation  Theory,”  IEEE  Trans,  informat.  Th.,  vol. 35,  pp. 1184-1205,  Nov. 
1989. 

[10]  S.  Y.  Rung  and  D.  W.  Lin,  “Optimal  Hankel  Norm  Model  Reductions:  Multi- 
variable  Systems,”  Int.  J.  Control,  vol. 39,  no. 6,  pp. 1115-1193,  1984. 

[11]  V.  M.  Adamjan,  D.  Z.  Arov,  and  M.  G.  Krein,  “Analysis  Properties  of  Schmidt 
Pairs  for  a  Hankel  Operator  and  the  Generalized  Schur-Takagi  Problem,”  Math. 
USSR.  Sbornik,  vol.  15,  no.l,  pp31-73,  1971. 

[12]  T.  Zhou  and  H.  Kimura,  “Time  domain  identification  for  robust  control,”  System 
&  Control  Letters,  vol. 20,  ppl67-178,  1993. 

[13]  A.  M.  King,  U.  B.  Desai  and  R.  E.  Skelton,  “A  Generalized  Approach  to  q- 
Markov  Covariance  Equivalent  Realizations  for  Discrete  Systems,”  Automatica, 
vol.24,  No.4,  pp.507-515,  1988. 

[14]  H.  Kimura,  “Positive  partial  realization  of  covariance  sequences,”  Modeling, 
Identification  and  Robust  Control,  C.  1.  Byrnes  and  A.  Lindquist.  Eds.  Ams¬ 
terdam:  North-Holland,  1987,  pp.  499-513. 


32 


[15]  P.  Dienes,  The  Taylor  Series,  New  York:  Dover  Publications,  1957. 

[16]  M.  Marden,  The  Geometry  of  the  Zeros  of  a  Polynomial  in  a  Complex  Variable, 
Amer.  Math.  Soc.,  Math.  Survey,  no.  3,  1949. 

[17]  S.  U.  Pillai,  T.  I.  Shim,  and  D.  C.  Youla,  “A  new  technique  for  ARMA- 
system  identification  and  rational  approximation,”  IEEE  Trans.  Signal  Process¬ 
ing,  vol.  SP-41,  no.  3,  pp.  1281-1304,  March  1993. 

[18]  W.  C.  Lee,  “Identification  of  Nonminimum  Phase  Stable  Systems,”  Ph.D  Disser¬ 
tation,  Department  of  Electrical  Engineering,  Polytechnic  University,  Brooklyn, 
New  York,  Jan  1995. 

[19]  D.  C.  Youla,  “The  EEE:  A  new  tunable  high-resolution  spectral  estimator,” 
Part  I,  Technical  note,  no.  3,  Department  of  Electrical  Engineering,  Polytechnic 
Institute  of  New  York,  Brooklyn,  New  York,  1980:  also  RADC  Rep.  RADC-TR- 
81-397,  AD  A114996,  February  1982. 

[20]  D.  C.  Youla,  “Two  observations  regarding  First- quadrant  Causal  BIBO  stable 
Digital  Filters,”  Proc.  IEEE,  vol.  78,  no.  4,  pp.  598-603,  April  1990. 


33 


(c)  Test  criteria  for  model  order  selection 


Figure  1:  Reconstruction  of  a  nonminimum  phase  ARMA(7,6)  model  from  its  partial 
impulse  response  sequence.  The  original  model  corresponds  to 

0.602  -  2.206^  +  4.160^2  -  4.952^^  +  4.3662:'*  -  2.7212:5  +  O.8882® 

~  1.0  -  3.5112  +  6.43822  -  8.05223  ^  7.87524  _  6.01825  +  3.1862^  -  0.8882^  ‘ 

Poles  of  H{z)  are  at  1.02Z0“ ,  l.OlZ  ±  35°  ,  l.OlZ  ±  51°  and  1.03Z100°. 

Zeros  of  H{z)  are  at  1.02Z  ±  22°  ,  0.769Z  ±  43°  and  1.05Z  ±  90° . 

The  reconstructed  model  is  given  by 

0.602  -  2.2O62  +  4.I6O22  -  4.95223  +  4.36624  -  2.72125  ^  0.8882® 

“  1.0  -  3.5112  +  6.43822  -  8.05223  7.875^4  _  6.01825  +  3.1862^  -  0.8882^  ’ 
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(a)  Original/reconst,  magnitude.  (b)  Original/reconst,  phase. 


Figure  2:  Reconstruction  of  an  MA(5)  model  from  its  partial  impulse  response  se¬ 
quence.  The  original  model  corresponds  to 

H{z)  =  0.8704  -  2.4714^  +  3.57822^  -  3.6560^^  -I-  2.0357^^  -  1.0^® . 

Zeros  of  H{z)  are  at  0.667Z0°  ,  0.952Z  ±  60°  and  1.2Z  ±  80°. 

The  reconstructed  model  is  given  by 

^  ,  0.8704  -  2.4714Z  +  3.57822^  _  3.6560^^  -h  2.03572^  -  l.Oz^ 

“  1  -f  2.1  X  10-1^2  -  1.5  X  10-1®22  -t-  4.1  X  10-1523  -b  2.6  X  10-14^4  +  7.2  X  10-1^25  ■ 
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dB 


TT 


0  7r/4  7r/2  37r/4 

(a)  Original/approx,  magnitude. 


0  7r/4  7r/2  37r/4  tt 

(b)  Original/approx,  phase. 


Figure  3:  Rational  approximation  of  a  nonrational  system.  The  original  nonrational 
transfer  function  is  given  by 

H{z)  =  H^{z)e-^  ^ 


where 

H  (z)  =  +  ^  u  =  _ 1.04 -1.312: +  2:^ _ 

’  1.02  -  1.75^  +  22’  >  1.09  -0.912  +  1.0622  +  1.0323  -0.862'!  +  25' 

The  stable  Fade  approximation  ARMA(11,11)  model  is  given  by  Hr{z)  =  Bm{z)/An{z), 
where 

An{z)  =  1  -  2.1112  +  2.35122  -  0.2652^  -  1.8692^  +  2.4592'^  -  1.0782® 

+0.1192^  +  0.2282®  +  0.05572^  +  0.00592^°  +  0.000262“  , 

Bm{z)  =  3.571  -  9.6092  +  15.72022  -  12.9742®  +  7.4402'‘  -  2.2432^  +  1.2272® 
-0.5552^  +  0.1262®  -  0.01882®  +  0.00162^°  -  0.000062“  . 

ZerosofA„(2)areatl.01Z±30°,  1.007Z±78.9° ,  1.0256Z±60°,  7.29Z±162°  ,  8.066Z+ 
45°  and  1.019Z180°  .  Zeros  of  Bm{z]  are  0.782Z  ±  53.2°,  1.107Z  ±  55°,  2.309Z  ± 
123°  ,  6.512Z  ±  67.5°  ,  10.17Z  ±  24°  and  3.37Z0.0°  . 
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(a)  Original/approx,  magnitude.  (b)  Original/approx,  phase. 


1,0  6,2  8,6  10,9  12,11  15,9  17,1  18,18  20,20 


(c)  Test  criteria  for  model  order  selection 

Figure  4:  Rational  approximation  of  a  nonrational  system  with  a  logarithmic  singu¬ 
larity  at  2  =  1.  The  original  nonrational  transfer  function  is  given  by 

0.797  -  1.036^  +  0.829^2 

“  1  -  1.4702  +  0.946z2  +  0.84123  _  1.27724  +  o.829^5  ~ 

Poles  of  E{z)  are  at  1.03Z  ±  39.6“ ,  1.04Z  ±  61.2“  and  1.05Z180°  .  Zeros  of  E{s) 
are  at  0.98Z  ±  50.4°  . 

The  stable  Fade  ARMA(14,4)  model  is  given  by  Er[z)  —  Bm{z)l An{z),  where 

Ar,{z)  =  1  - -2.85362  +  3.3370^2-  0.95902^ -2.1421.2^+2.9183;2®- 1.5714/ 
+0.2655220.01882* +  0.00532^0.00222^°0.0011z“  +5.21  x  10-‘‘222 
+2.11  X  10"^2^*  +  5.35  X  10-^22'* , 

Bm{z)  =  -0.0009  -  0.79552  +  1.739822 -  1.74552* +0.733052/ 

Zeros  of  Bm{z)  are  at  0.98Z  ±  50.4° ,  1.1323Z0.0°  and  0.0011Z180.0° . 
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3.  Stacked  Piei^ocermaic  Actuators  for  Active  Damping 


3.1  Introduction 

Utilization  of  active  materials  for  vibration  damping  is  recently  receiving  increased 
attention  [1-6,25].  Research  and  development  of  completely  integrated  structures 
with  embedded  actuators,  sensors,  signal  processing,  and  control  systems  are  of  prime 
interest.  Different  electroactive  materials  such  as  piezoceramics,  shape  memory  alloys, 
elect rorheologi cal  fluids,  polymer  biomaterials,  and  magnetorestrictive  materials  are 
being  considered  for  such  purposes.  Piezoelectric  ceramics  are  effective  distributed 
strain  actuators  as  well  as  sensors  due  to  their  high  stiffness,  good  linearity,  ease  of  in¬ 
tegration,  low  temperature  sensitivity,  and  relatively  low  noise.  In  this  report,  stacked 
piezoceramics  are  mainly  considered  for  actuation  purposes.  Also,  piezoceramics  with 
interdigitated  electrodes  are  considered  [8,9].  Dynamical  models  of  stacked  piezoce¬ 
ramics  are  reported.  Several  experimental  setups  (ranging  from  beams  and  plates  to 
multi-link  flexible  manipulators)  have  been  developed  to  study  modeling  and  control 
of  smart  structures  at  the  CRRL. 

One  of  the  main  issues  in  smart  structures  is  the  control  system.  Foremost  in  the 
design  of  advanced  control  systems  is  the  need  to  satisfy  performance  specifications. 
In  fact,  one  often  sacrifices  stability  margin  or  accepts  conditionally  stable  systems  in 
order  to  realize  superior  performance.  So  one  should  begin  with  an  optimal  solution 
according  to  some  measure.  In  the  present  work,  an  H2  measure  is  considered.  When 
the  optimal  solution  exceeds  specifications,  the  cost  functional  can  be  increased  and 
the  increment  traded  for  an  improvement  in  stability  margin.  That  is,  one  should 
maximize  stability  margin  subject  to  a  constraint  on  the  H2-norm  cost  increment. 
When  the  stability  margin  obtained  in  this  way  is  not  adequate,  then  one  must  use  a 
different  controller  for  different  operating  profiles  of  the  plant:  typically,  one  should 
not  give  up  performance  just  for  the  sake  of  operating  with  a  single  fixed  controller 
when  high  performance  is  required.  In  order  to  attain  an  analytical  solution  to  this 
fundamental  tradeoff  problem,  the  approximate  measure  of  stability  margin  chosen  is 
also  a  quadratic  cost  function.  A  complete  derivation  of  the  results  for  the  advocated 
control  system  is  given  in  [20]. 

This  new  control  design  methodology  for  tradeoffs  between  stability  margin  and 
performance  is  applied  to  vibration  damping  and  pointing  of  smart  structures.  These 
systems  do  provide  nice  benchmark  problems  in  control  system  design.  The  afore¬ 
mentioned  control  design  methodology  is  applied  to  a  flexible  pointing  system  with 
embedded  active  materials  and  experimental  results  are  presented  to  show  the  efficacy 
of  the  control  design  methodology  and  the  added  advantage  of  utilization  of  active 
materials. 

3.2  Stacked  Piezoceramic  Elements 

In  this  section,  the  dynamical  model  of  stacked  piezoceramic  elements  (e.g.,  plates, 
discs,  etc.)  attached  electrically  in  parallel  is  developed.  A  typical  stack  of  discs  is 


.38 


F 


Figure  3.1:  A  stack  of  N-piezoceramic  discs. 


shown  in  Figure  3.1  where  the  piezoceramic  elements  are  assumed  to  be  polarized 
through  their  thickness.  The  electric  work  done  by  each  element  of  the  stack  is 
derived  to  generate  the  dynamical  model  and  to  calculate  the  force  generated  by  the 
stack. 

The  strain  through  the  piezoelectric  ceramic  is  expressed  as 
S=[Si  S2  Sz  S4  Ss  Se]'  =  L^u  (3.1) 


where  (•)'  denotes  the  transpose,  is  a  6x3  differential  operator,  and  u  —  [u^  Uy 
is  a  3  X  1  displacement  vector. 

The  constitutive  relation  between  the  stress  (T),  strain  (S),  electric  field  (E),  and 
electric  charge  distribution  (D)  in  each  of  the  piezoceramic  elements  is  given  by 


T 

D 


(3.2) 


where  is  the  stiffness  matrix  of  the  piezoceramic  taken  at  constant  electric  field, 
is  the  matrix  of  dielectric  constants  taken  at  constant  strain,  and  e  is  the  matrix 
of  material  constants  of  the  piezoceramic.  All  quantities  are  taken  with  respect  to  a 
local  coordinate  system  pertaining  to  the  direction  of  polarization. 

The  kinetic  and  potential  energies  of  an  element  in  the  piezoceramic  stack  are 
given  respectively  as 

T  —  ~~ps  f  iiii  df],  and  {/  =  -  /  S’T 

j-i  J ^ g  ^ 

where  is  the  mass  volume  density  of  the  piezoceramic  element  and  ilg  denotes  the 
spatial  boundary  of  the  section  of  the  piezoceramic. 


(3.3) 
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Replacing  T  by  its  equivalent  from  Eq.(3.2)  simplifies  the  expression  for  the  po¬ 
tential  energy  into 

[7  =  1/  (S'c^S  -  S'e'E)d^.  (3.4) 

2  JUs 

The  electric  work  done  by  the  charge  in  each  element  of  the  stack  is  expressed  as 
We  =  \-  j  E'D  dn.  (3.5) 

2  JQ^ 


Replacing  D  in  the  above  equation  by  its  value  in  Eq.  (3.2)  yields 

We  =  -  [  E'eS  dn  +  -  f  E'e^E  d^  (3.6) 

where  the  first  term  signifies  the  actuation  capability  and  can  be  expanded  as 

^act  =  ^  /  (esiS'i  +  632^2  +  633^3)  dO.  (3.7) 

iuTjfp  •/O3 

where  the  first  subscript  refers  to  the  direction  of  the  electric  field  and  the  second 
refers  to  the  direction  of  the  strain.  The  3-direction  is  the  direction  of  polarization 
(here  corresponding  to  the  x-direction),  the  2-direction  corresponds  to  z  and  the  1- 
direction  corresponds  to  y.  Note  that  54  =  S's  =  =  0  signifying  the  absence  of 

shear  stresses. 

For  brevity,  the  piezoceramic  elements  are  assumed  to  be  identical  and  of  a  simple 
geometry  (e.g.,  discs  or  plates)  in  the  following  analysis.  Furthermore,  the  average 
electric  field  through  each  of  the  elements  can  be  expressed  as' 


E  - 


0  0  f 

tp 


(3.8) 


where  V  is  the  applied  voltage  and  tp  is  the  thickness  of  the  piezoceramic  element. 

The  dominant  strain  is  in  the  direction  of  polarization  (i.e.,  the  3-  direction). 
Hence,  one  can  assume  that  only  Ss  is  non-zero.  It  then  follows  from  Eqs.  (3.4)  and 
(3.8)  that 


U  = 


1 


A 

/i.pc.33 


ra+tp 

J  a 


S3  dx 


VApCsi  /■“+'? 


2tr 


f 


S3  dx 


“  2  ””7  iax 


Ua;(x,i)]  dx 


VApe3i  f°-+h  7  d 
I  dx 


A-pCsi  r 

2tp  Ja 


Ux{x, t)  )  dx 


(3.9) 


^We  consider  uniform  field  through  the  piezoceramic  since  thin  elements  are  normally  utilized. 
However,  the  non-uniform  field  may  similarly  be  handled. 
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where  a  and  a  +  are  the  spatial  limits  of  the  piezoceramic  element  along  the  x-axis 
and  Ap  is  the  area  of  the  faces  of  each  of  the  piezoceramic  elements  constituting  the 
stack. 

Similarly,  the  component  of  displacement  of  largest  significance  is  Us  (x-direction); 
therefore,  by  neglecting  the  displacements  in  other  directions,  the  kinetic  energy  sim¬ 
plifies  into 

T  =  -psAp  ul{xA)  dx  (3.10) 

Finally,  the  total  kinetic  and  potential  energies  of  the  stack  are  obtained  by  summing 
up  those  of  the  individual  piezoceramic  elements  given  in  Eqs.  (3.9)  and  (3.10).  There¬ 
after,  Hamilton’s  principle  may  be  invoked  to  derive  the  dynamics  of  the  piezoceramic 
stacks^.  Furthermore,  a  discretization  of  the  distributed  parameter  model  may  be  per¬ 
formed  to  obtain  a  lumped  parameter  description  of  the  piezoceramic  stacks.  This 
may  be  achieved  through  an  assumed  mode  approach  on  the  displacement,  i.e., 

m 

Ux{x,t)  (3-11) 

t  =  l 


where  m  is  the  number  of  retained  modes,  the  d^i{x)  are  the  mode  shape  functions 
and  the  rii(t)  are  the  amplitude  functions. 

The  following  analysis  establishes  the  level  of  force  such  piezoceramic  stacks  may 
generate.  Neglecting  all  strains  in  other  directions  than  the  polarization  direction, 
the  work  done  by  the  actuator  given  in  (3.7)  can  be  simplified  to 

=  r'^S^dx.  (3.12) 

Ztp  J  a 

Approximating  the  strain  through  the  length  as  a  constant,  it  can  be  expressed  as 

=  ^  (3.13) 

where  L  is  the  total  length  of  the  stack  and  AL  is  the  change  in  length  of  the  stack 
due  to  an  applied  electric  field.  Then  the  expression  for  the  actuation  capability  of 
each  element  becomes 


IT. 


act 


V  ApC^sAL 
21 


Remark:  the  total  actuation  of  the  stack  is  N  x  Wact- 


(3.14) 


^Tip  mass  energy  and  spring  energy  due  to  preload  may  also  be  added. 
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To  compute  the  force  generated  by  one  element  of  the  stack,  consider  the  mechan¬ 
ical  work  done: 


II 

(3.15) 

where  F  —  [Fx  Fy  F^]  is  a  three  dimensional  vector  with 

Fx  =  [  f  Tzdydz, 

JQz 

(3.16) 

Fy  —  f  f  Tidxdz, 

JQz  JQx 

(3.17) 

Fx  =  [  f  T2dxdy, 

(3.18) 

and  and  are  the  spatial  domains  of  the  piezoceramic  element  along  the  x. 

y  and  z  axes  respectively.  In  Eq.  (3.15),  u  is  the  displacement  vector  which  can  be 
expressed  in  terms  of  the  normal  strains  (in  its  local  coordinates) 


r 

= 

Szx  S\y  S2Z 


f 


(3.19) 


When  the  strain  components  are  constant  along  their  respective  directions  throughout 
the  piezoceramic  stack,  du  can  be  expressed  as 


du  =  [  Szdx  Sidy  S2dz  ]'•  (3.20) 

Thus,  the  mechanical  work  becomes 


W 


—  f  FxSsdx  +  I  FyS\dy  +  /  FzS2dz 

J  (X  J  -y  J ^ 


(3.21) 


Since  the  component  of  the  force  being  considered  (x-component)  is  much  larger  than 
the  other  two  force  components,  Eq.  (3.21)  may  be  approximated  by 

ra-\-tp 

W=  F^Szdx-  (3.22) 

J  a 

Having  a  constant  force  along  the  x-direction,  it  may  be  computed  by  equating 
the  actuation  energy  to  the  mechanical  energy  to  obtain 

T’x  =  \A,Eezz  =  (3.23) 


This  is  the  force  due  to  each  of  the  elements  in  the  stack;  therefore,  the  total  force 
exerted  by  the  stack  is 


F.= 


NVA.ezz 

2tp 


(3.24) 
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Figure  3.2:  A  piezoceramic  plate  with  interdigitated  electrodes. 


If  the  stack  is  preloaded  with  a  spring  of  stiffness  c^,  633  in  Eq.  (3.24)  changes  to  633 


where 

.  ct 

633  —  ,  633 

Ct  +  Cs 

(3.25) 

and 

NAp 

Ct  =  C33  ^  • 

(3.26) 

3.3  Piezoceramic  Plates  With  Interdigitated  Electrodes 

In  this  section,  a  model  of  a  piezoceramic  plate  with  interdigitated  electrodes  is 
presented.  Fig.  3.2  gives  the  top-view  of  a  piezoceramic  plate  with  interdigitated 
electrodes.  The  interdigitated  piezoceramic  is  a  special  case  of  the  aforementioned 
stacked  actuators.  The  energy  expressions  are  needed  to  compare  the  actuation  capa¬ 
bilities  of  the  piezoceramic  with  interdigitated  electrodes  to  those  of  a  piezoceramic 
with  face-plane  (conventional)  electrodes. 

Taking  advantage  of  the  symmetry  of  the  interdigitated  pattern,  one  first  considers 
a  section  of  the  piezoceramic  including  one  electrode  of  each  polarity  (Figure  3.3.) 
Furthermore,  since  the  thickness  of  the  section  is  small  with  respect  to  the  spacing 
between  the  electrodes,  the  transverse  component  of  the  electric  field  is  neglected 
and  an  average  value  of  the  longitudinal  component  will  be  considered.  The  electric 
potential  at  all  points  under  an  electrode  is  equal  to  that  of  the  electrode;  therefore, 
no  electric  field  is  considered  in  these  areas.  Hence,  the  average  electric  field  through 
a  section  lies  in  the  area  between  the  electrodes  and  is  given  hy  E  ~  j  where  V  is  the 
applied  voltage  to  the  electrodes  and  g  is  the  distance  between  electrodes  of  opposite 
polarity. 
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Figure  3.3;  A  section  of  the  interdigitated-electrode  piezoceramic. 


The  assumptions  made  here  reduce  the  piezoceramic  with  interdigitated  electrodes 
to  a  stack  of  identical  piezoceramic  sections  with  each  section  as  shown  in  Fig  3.3. 
Therefore,  the  results  of  Section  3.2  are  applicable  here. 

Since  thin  piezoceramic  sections  are  of  main  interest,  the  dominant  strain  is  in 
the  3-direction  (along  the  length);  therefore,  from  Eq.  (3.12)  the  actuation  energy  of 
each  section  is  given  by 


IF..*  = 


Vwststzz  /'“+® 


(3.27) 


where  and  4  are  the  width  and  thickness  of  the  piezoceramic  section  and  a  and 
a  +  g  are  the  spatial  limits  of  the  section  along  the  x-axis.  The  expression  for  the 
potential  energy,  kinetic  energy,  and  force  generated  by  JV  piezoceramic  sections  is  as 
given  in  Eqs.  (3.9),  (3.10),  and  (3.24),  respectively,  with  Ap  replaced  by  Wgig- 

To  compare  the  actuation  level  of  a  piezoceramic  with  interdigitated  electrodes 
to  those  of  a  conventional  piezoceramic  plate  with  face  plane  electrodes,  consider 
two  piezoceramic  plates  with  the  same  length,  width  and  thickness.  Ip,  Wp,  and  tp, 
respectively.  In  the  case  of  the  interdigitated  piezoceramic,  the  dominant  strain  is 
in  the  3-direction  whereas  the  dominant  strain  of  the  conventional  piezoceramic  is  in 
the  1-direction.  Comparing  the  actuation  of  the  proposed  piezoceramics  yields 


IFge*. 


V Wsts  Y" 

Vwpt 

2tp 


633^3 

^  631 


(3.28) 


where  Wacti  is  the  actuation  level  of  the  ceramic  with  interdigitated  electrodes,  Wactc 
is  the  capability  of  the  conventional  piezoceramic,  N  is  the  number  of  section  in  the 
interdigitated  piezoceramic,  and  ai  is  the  spatial  limit  of  the  piezoceramic  section. 
Normalizing  electric  fields  and  taking  tp  —  2ts  simplifies  the  expression  into 


IFgc*.  _W,e3sE^lfaT^S3 

Wact,  2wpe3:f:+‘^Si 


(3.29) 


Similarly,  comparing  the  forces  resulting  from  both  ceramics  yields 
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Ratio  of  Actuation  Enorgias  with  Normalized  Electric  Field 


Figure  3.4:  Comparison  of  actuation  levels  due  to  conventional  and  interdigitated 
piezoceramics. 


Fi  _  NwstsEczz 


(3.30) 


where  Fi  is  the  force  exerted  by  the  interdigitated  piezoceramic  and  Fc  is  the  force 
due  to  the  conventional  piezoceramic.  Note  that  for  a  normalized  field  and  tp  =  2ts, 
the  ratio  simplifies  to 


a  ^  Nw,e,s  p 

F^  2uip63i 

To  compare  the  performance  in  vibration  damping  of  the  two  proposed  piezo¬ 
ceramic  plates,  attach  each  of  them  to  an  identical  clamped-free  beam  at  identical 
positions.  Fix  the  width  tg  of  the  fingers  of  the  interdigitated  electrodes  and  let 
tg  =  tp.  The  actuation  energies  of  the  two  piezoceramic  plates  are  compared  us¬ 
ing  Eq.  (3.29)  and  the  mode  shape  functions  of  the  clamped-free  beam  (normalized 
electric  field).  Figure  3.4  shows  the  ratio  Wacti/Wactc  for  different  values  of  g/tg. 


3.4  Wiener-Hopf  Design 

An  analytical  solution  to  the  tradeoffs  of  performance  versus  stability  margin  has 
been  attained  and  reported  in  [20]  (in  these  proceedings).  The  approach  is  based  on 
H2  measures  of  performance  and  stability  margin.  Furthermore,  a  novel  approach 
is  also  introduced  which  allows  structured  perturbations  in  the  coprime  polynomial 
matrix  fraction  description  of  the  plant  transfer  matrix  to  be  taken  into  account.  The 
reader  is  referred  to  [20]  for  details. 
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Simulation  ol  Optn  and  Cbsed-Loofi  Raaponsai 


Figure  3.5:  Simulation  results  for  open  and  closed-loop  system. 


3.5  Simulation  and  Experimental  Results 


An  interdigitated  electrode  pattern  was  deposited  onto  each  side  of  a  PZT,  Pb(Zr,  TI)03, 
plate.  The  PZT  plate  was  of  type  5H,  manufactured  by  Morgan  Matroc,  Inc.  The 
electrodes  were  deposited  onto  the  PZT  by  evaporation.  First,  the  area  which  was 
to  remain  unelectroded  was  covered  with  a  stainless  steel  mask.  We  assured  that  the 
electrodes  on  both  sides  of  the  PZT  plate  were  aligned  through  a  registration  scheme. 
Thereafter,  a  layer  of  nickel  was  deposited  on  each  side  of  the  ceramic;  nickel  provided 
a  strong  bond  to  the  ceramic  whereas  the  conductivity  of  the  resulting  electrode  was 
not  sufficient.  To  improve  the  conductivity  of  the  electrodes,  a  layer  of  gold  was 
deposited  over  the  nickel  resulting  in  a  maximum  resistance  of  15  Ohms  along  any 
pair  of  electrodes.  Electrical  contact  was  established  through  the  use  of  a  two  part 
conductive  glue. 

A  controller  based  on  the  Wiener-Hopf  {-H2)  design  methodology  is  designed 
for  a  clamped-free  beam  with  an  interdigitated  piezoceramic  actuator  and  a  con¬ 
ventional  piezoceramic  sensor.  Perfect  bonding  between  the  piezoceramic  actuator 
(2.5"  X  1.5"  X  0.020")  and  the  rectangular  aluminum  beam  (17"  x  5"  x  1/16")  is  as¬ 
sumed,  so  that  the  overall  model  is  not  affected  by  the  adhesive  layer.  The  sensor  is 
mounted  further  above  the  actuator  and  due  to  its  non-collocation  with  the  actuator, 
a  non-minimum  phase  plant  is  attained.  The  analytical  modeling  and  the  experi¬ 
mental  results  do  match  with  each  other.  The  parameters  in  the  Wiener-Hopf  design 
were  chosen  such  that  the  advocated  controller  was  not  proper.  The  controller  has 
four  zeros  and  three  poles.  This  was  done  in  order  to  achieve  velocity  feedback  for  a 
better  active  damping  of  the  structure.  The  controller  is  as  follows: 


0.188s^  -h  63.56s2  -f-  3754. +  387420s  4-  3560400 
s3  -I- 171. Is2  -1-  4997.7s  -|-  11452 


(3.32) 


The  open  and  closed-loop  responses  of  the  structure  for  an  impulse  type  force  distur¬ 
bance  at  the  tip  of  the  beam  is  shown  in  Fig.  3.5.  The  applied  electric  field  does  not 
exceed  15  Kvolts/in. 
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Expenmenl*)  Oosad  and  Open-Loop  Responses 


Figure  3.6:  Experimental  open  and  closed-loop  responses. 

The  obtained  controller  was  tested  on  the  experimental  setup  at  the  CRRL.  The 
controller  was  discretized  and  implemented  by  a  ’486  based  PC  with  a  sampling 
period  of  2  msec.  The  DAC  output  was  filtered  with  a  second-order  filter  having  a 
natural  frequency  of  150  Hz.  A  high  voltage  amplifie  was  used  to  step  up  the  control 
voltage  to  the  piezoceramic  actuator.  The  open  and  closed-loop  responses  are  shown 
in  Fig. 3. 6.  A  significant  improvement  in  the  settling  time  of  the  system  has  been 
achieved.  We  are  currently  implementing  controllers  on  a  flexible  pointing  system. 

3.6  Conclusions 

Dynamical  models  of  stacked  piezoceramic  actuators  were  established.  The  stan¬ 
dard  piezoceramic  stack  actuators  and  the  interdigitated  electrode  patterns  were  con¬ 
sidered.  The  stacked  piezoceramics  were  utilized  for  active  structural  damping.  A 
new  robust  control  design  methodology  based  on  H2  design  was  advocated  and  shown 
to  be  an  effective  design  methodology. 
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